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Abstract 

We analyze the global dynamics of Bianchi type I solutions of the Einstein equations with 
anisotropic matter. The matter model is not specified explicitly but only through a set of mild 
and physically motivated assumptions; thereby our analysis covers matter models as different 
from each other as, e.g., coUisionless matter, elastic matter and magnetic fields. The main 
result we prove is the existence of an 'anisotropy classification' for the asymptotic behaviour 
of Bianchi type 1 cosmologies. The type of asymptotic behaviour of generic solutions is de- 
termined by one single parameter that describes certain properties of the anisotropic matter 
model under extreme conditions. The anisotropy classification comprises the following types. 
The convergent type A_|_: Each solution converges to a Kasner solution as the singularity is 
approached and each Kasner solution is a possible past asymptotic state. The convergent types 
B_|_ and C-I-: Each solution converges to a Kasner solution as the singularity is approached; 
however, the set of Kasner solutions that are possible past asymptotic states is restricted. 
The oscillatory type D-|-: Each solution oscillates between different Kasner solutions as the 
singularity is approached. Furthermore, we investigate non-generic asymptotic behaviour and 
the future asymptotic behaviour of solutions. 



1 Introduction 

A pivotal feature in the study of spatially homogeneous cosmologies is the fact that the (asymp- 
totic) dynamics of cosmological solutions of the Einstein equations is largely determined by the 
spatial geometry of the initial spacelike hypersurface. This is reflected in the systematics of the 
Bianchi classification of spatially homogeneous cosmologies. In contrast, a systematic analysis of 
the influence of the matter model on the (asymptotic) dynamics of solutions is not available at 
present. It is known that the dynamics of solutions (and the asymptotics towards the initial singu- 
larity in particular) strongly depends on the matter source that is considered; the relatively simple 
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behaviour of the vacmini and the perfect fluid case, see, e.g., [22]. is replaced by a considerably 
more intricate behaviour in the case of anisotropic matter models [3, 6, 11, 15, 19]. The purpose 
of this paper is to analyze systematically and in detail in which way anisotropies of the matter 
source determine the dynamics of cosmological solutions. 

In our analysis we consider spatially homogeneous solutions of Bianchi type I. The reason for 
this choice is that Bianchi type I is the 'foundation' of the Bianchi classification. Generally, 
the dynamics of the 'higher' types in the Bianchi classification is based on the dynamics of the 
'lower' types; for instance, the asymptotics of solutions of Bianchi types VIII and IX can only 
be understood in terms of solutions of Bianchi type I and II, see [7] for an up-to-date discussion. 
Likewise, our understanding of cosmological models in general is conjectured to be built on the 
dynamics of spatially homogeneous models and thus on the dynamics of Bianchi type I models in 
particular, see [8] and references therein. 

The dynamics of Bianchi type I cosmological models where the matter is a perfect fluid is well- 
known [22]. In this context, it is customary to assume that the the perfect fluid is represented 
by an energy density p and pressure p that obey a linear equation of state p = wp with w = 
const e (—1,1). Each Bianchi type I perfect fluid solution isotropizes for late times; this means 
that each solution approaches a (flat) Priedmann-Robertson- Walker (FRW) solution for late times. 
Towards the initial singularity, each solution approaches a Kasner solution (Bianchi type I vacuum 
solution), i.e., to leading order in the limit t ^ Q the metric is described by 

ds^ = -dt^ + ai t'^P'dx^ (g) dx^ + t'^^^dx'^ (g) dx'^ + as t'^^'^dx^ (g) dx^ , 

where ai, a2, are positive constants, andpi, p2, ps are the Kasner exponents, which are constants 
that satisfy P1+P2+P3 = Pi +P2 +P3 ~ 1- Conversely (and importantly), for each Kasner solution 
(including the Taub solution, for which {pi,p2,p:i) = (1,0,0) and permutations) there exists a 
Bianchi type I perfect fluid model that converges to this solution as t ^ 0; in other words, each 
Kasner solution is a possible past asymptotic state. 

In this paper we are concerned with Bianchi type I cosmological models with anisotropic matter. 
We do not specify the matter model explicitly; on the contrary, the assumptions that we impose on 
the matter source are so mild that we can treat several matter models at the same time, examples 
being matter models as different from each other as coUisionless matter, elastic matter (for a wide 
variety of constitutive equations) and magnetic fields (aligned along one of the axes). The main 
results we derive are presented as Theorems 1-3 in Section 6 and summarized in Figure 7, but let 
us give a rough non-technical overview: 

First consider the asymptotic behaviour of models towards the future. For 'conventional' matter 
sources — where 'conventional' means that the isotropic state of these anisotropic matter models is 
energetically favorable — we find that each associated Bianchi type I solution isotropizes towards 
the future, i.e., each solution approaches a FRW solution for late times. However, there also exist 
matter models whose isotropic state is unstable. In this case, the future asymptotic behaviour 
is completely different; isotropization occurs but it is non-generic; generically, solutions approach 
(lifFcrent self-similar solutions, e.g., there exist solutions that approach Kasner solutions as f — > 00. 
It should be emphasized that these models are often compatible with the energy conditions. 

Second consider the asymptotic behaviour towards the singularity. We observe a rather diverse 
asymptotic behaviour, but the past dynamics of generic Bianchi type I solutions is again inti- 
mately connected with the Kasner states. Interestingly enough, the details of the past asymptotic 
dynamics of solutions are governed by one particular parameter, (3, that describes certain 'asymp- 
totic properties' of the matter (the properties of the matter under extreme conditions — extreme 
stress, etc. — that are found close to the spacetime singularity); this parameter can also be regarded 
as a measure for the degree of (dominant) energy condition violation under extreme conditions. 
(Clearly, we always assume that the dominant energy condition is satisfied under 'normal condi- 
tions' of the matter.) If the matter satisfies the dominant energy condition also under extreme 
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conditions, then each Bianchi type I solution approaches a Kasner sohition, and conversely, each 
Kasner solution is a possible past asymptotic state; the Kasner solutions are on an equal footing in 
this respect. (The case where the energy condition is satisfied only marginally is slightly different; 
the Taub solutions play a special role in that case.) If the dominant energy condition is violated 
under extreme conditions, this result breaks down; however, we must distinguish between 'under- 
critical' violations of the energy condition (corresponding to a small value of /3) and 'over-critical' 
violations (corresponding to a large value of /?). For an under-critical violation of the energy con- 
dition, each Bianchi type I solution approaches a Kasner state, but the converse is false, i.e., there 
exist Kasner solutions that are excluded as possible past asymptotic states. Only Kasner solutions 
that are 'sufficiently different' from the Taub solution, i.e., only Kasner solutions with {pi,P2,P3) 
sufficiently different from (1,0,0) and permutations, are attractors for Bianchi type I solutions. 
The set of potential asymptotic states is smaller if the degree of the energy condition violation 
is larger. Finally, for an over-critical violation of the energy condition, Bianchi type I solutions 
with anisotropic matter do not converge to a Kasner state as t ^ 0, but they follow a sequence of 
Kasner states ('epochs'), in each of which their behaviour is approximately described by a Kasner 
solution. This is a type of behaviour that resembles the Mixmaster oscillatory behaviour that is 
expected for Bianchi type VIII and IX (vacuum) cosmologies (and generic cosmological models); 
note, however, that both the details and the origin of these oscillations are very different — in the 
present case, the oscillatory behaviour is not related to the Mixmaster map; it is a consequence of 
(over-critical) energy condition violation; in the Mixmaster case, the reason is the geometry of the 
problem. Independently of the energy conditions, there also exist a variety of Bianchi type I solu- 
tions with anisotropic matter that exhibit a past asymptotic behaviour that is not connected with 
Kasner asymptotics; however, these are non-generic solutions; we prove that there exist essentially 
three different types of non-generic asymptotics. 

The outline of the paper is as follows. In Section 2 we discuss the Einstein equations in Bianchi 

type I symmetry and introduce our main assumption on the anisotropic matter model: The stress- 
energy tensor is assumed to be represented by a function of the metric and the initial data of the 
matter source. We shall restrict ourselves to the study of diagonal models, i.e., solutions for which 
the metric is diagonal in the standard coordinates of Bianchi type I symmetric spacctimes. In 
Section 3 we rewrite the Einstein equations in terms of new variables and show that the essential 
dynamics is described by a system of autonomous differential equations on a four dimensional 
compact state space — the reduced dynamical system. Section 4 is by far the longest and most 
technical part of the paper. It contains a detailed analysis of the flow induced by the reduced 
dynamical system on the boundary of the state space. In Section 5 we study the local stability 
properties of the fixed points of the reduced dynamical system. Finally, in Section 6 we present 
the main results: Theorems 1-3 and Figure 7. The proofs are based on the results of Sections 4 
and 5 and make use of techniques from the theory of dynamical systems. The purpose of Section 7 
is to bring to life the general assumptions on the matter source that are made in Sections 2-5. We 
give three important examples of matter models to which our analysis applies straightforwardly: 
collisionless (Vlasov) matter, elastic matter and magnetic fields. Finally, Section 8 contains some 
concluding remarks together with an outlook of possible developments of the present work. 



2 The Einstein equations in Bianchi type I 

We consider a spatially homogeneous spacetime of Bianchi type I. The spacetime metric can be 
written as 

ds^ = -df + gij{t)dx'dx^ (ij = 1, 2, 3) , (1) 

where gij is the induced Riemannian metric on the spatially homogeneous surfaces t = const. By 
kij we denote the second fundamental form of the surfaces t = const. The energy-momentum tensor 
(m, = 0, 1,2,3) contains the energy density p = Too and the momentum density jk = Tofe. 
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Einstein's equations, in nnits c = 1 = SttG, decompose into evolution equations and constraints. 
The evolution equations are 

dtgij = -2guk^j , dtk'j = tr fc fc^ - T) + ^S'j (tr T - p) , (2a) 
the Hamiltonian constraint reads 

{tvkf -k'jk\-2p = 0, (2b) 

and the momentum constraint is jk = (which is due to the fact that Bianchi I spacetimes are 
spatially flat). The latter does not impose any restriction on the gravitational degrees of freedom, 
which are represented by gij and fc'^-, but only on the matter fields. In general, the system (2) 
must be complemented by equations describing the evolution of the matter fields. 

For certain matter models, the energy-momentum tensor takes the form of an explicit functional of 
the metric gij. Let us give some examples. For a (non-tilted) perfect fluid, the energy-momentum 
tensor reads Tf^i.dx^dx" = pdt^ + p gijdx'^dx^ , where the energy density p and the pressure p are 
connected via a barotropic equation of state p{p). When we define n as n = (detg)-i/2, the 
Einstein equations (2) imply p{p) + p — n{dp/dn) (which is a special case of equation (3) below). 
Hence p is obtained as a function p = p{n). Clearly, prescribing p as a function of n is equivalent to 
prescribing an equation of state for the fluid (modulo a scaling constant). Likewise, for coUisionless 
matter the components of the energy- momentum tensor are given as functions of the metric gij, 
which depend on the initial data for the matter (which is described by the distribution function of 
particles in phase space; see Section 7 for details). We thereby obtain the energy density p and the 
principal pressures Pi (which arc defined as the eigenvalues of T*^) as p = p{gij) and Pi = Pi(gjk)- 
Finally, for elastic matter the constitutive equation of state of the material determines p and Pi as 
functions of gij ; for details we refer again to Section 7. Motivated by these examples we make our 
fundamental assumption on the matter model. 

Assumption 1. The components of the energy- momentum tensor are represented by smooth (at 
least C^) fmictions of the metric gij. We assume that p = p{gij) is positive (as long as gij is 

non-degenerate) . 

Remark. We note that the particular form of the functions that represent the components of the 
energy-momentum tensor may depend on a number of external parameters or external functions 
that describe the properties of the matter, and/or the initial data of the matter field(s); see 

Section 7 and [21] for examples. 

By Assumption 1, the evolution equations of the matter fields in Bianchi type I are contained in 
the Einstein evolution equations (2a) via the contracted Bianchi identity. 

Remark. The regularity of T^j as a function of the spatial metric ensures that the Cauchy problem 
for the evolution equations with initial data at t = Iq is locally wcU-poscd on a time interval 
{t-,t+), t- < to < t+. The requirement p > and the constraint equation (2b) imply that the 
mean curvature tr k never vanishes in the interval of existence. Without loss of generality we assume 
that trfc < on (t_,f-|-), i.e., we consider an expanding spacctime; the case trfc > is obtained 
from tr A; < by replacing t with —t. If the energy momentum tensor satisfies the inequality 
trT -|- 3yO > (which corresponds to Assumption 3 in Section 3), then the mean curvature is non- 
decreasing, the solutions of (2) exist in an interval (i_,+oo) with i_ > — oo, and the spacctime 
is future geodesically complete. By a time translation we can assume t- = 0. Under additional 
conditions on the limit f — > corresponds to a curvature singularity. We refer to [16, 18] 
for a proof of these statements in the context of general spatially homogeneous cosmologies; for 
the special case we consider these results follow straightforwardly from our formulation of the 
equations, cf. (11). 
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Lemma 1. Assumption 1 implies 

In particular, the functional dependence of the energy-momentum tensor on the metric Qij is com- 
pletely determined by the function p{gij)- 



Proof. Differentiating (2b) w.r.t. t and using the evolution equation for fc'^- we obtain 

5tp = fc^,(TV+(5>). (4a) 
On the other hand, using p = p{gij) and the evolution equation for the metric we obtain 



Equating the r.h. sides of (4a) and (4b) gives an identity between two polynomials in k■'^ with 
coefficients that depend only on the metric gij; therefore the corresponding coefficients must be 
equal, which gives (3). □ 



Equation (3) will play a fundamental role in the following. We emphasize that the energy- 
momentum tensor miist be independent of the second fundamental form for the proof of Lemma 1 
to hold. We refer to Section 7 for a derivation of (3) within the Lagrangian formalism. 

Initial data for the Einstein-matter system are given by gij{to), k''j(to) and the initial values of the 
matter fields. Without loss of generality we can assume that gij{to) and k^j{to) are diagonal (by 
choosing coordinates adapted to an orthogonal basis of eigenvectors of k^j{to))- 

Assumption 2. We assume that there exists initial data for the matter such that T*^- is diagonal in 

any orthogonal frame. 

Uniqueness of solutions of the evolution equations then implies that (g-y , fc*^ , T*^ ) remain diagonal 
for all times. We refer to solutions of this type as diagonal models; henceforth we restrict ourselves 
to these models. 



3 Dynamical systems formulation 



We introduce a set of variables for which the Einstein evolution equations (2a) decouple. The 
reduced system we thereby obtain is a dynamical system on a bounded state space. Let 

X = g^^ + g^^ + g^^ (dimensional variables), 

gii 

Si = — (dimensionless variables). 

X 

There is no summation over the index i in these definitions. The division by H is not a restriction, 
since H > for all solutions, cf. the remark in Section 2. The dimensionless variables satisfy the 
constraints 

Si + 112 + S3 = , .Si + S2 + S3 = 1 . 

The transformation from the six variables {g",k^i) to the variables {H,x, Si,J^i) that satisfy 
the constraints is one-to-one. The variable x can be replaced by n = (det(?)~^/^, since x = 
n2/3(siS2S3)-^/'. 



H 



trfc 
"~3~ 

k\ 



5:. = ---l, 



5 



We define dimensionless matter quantities by 

p P,. T\ 



Wi = — 

p 



p 



1 / \ p 

W = -{Wi +W2+ W3) = - 

3 p 



tiiere is no summation over i. Tlie quantity p is th.e isotropic pressure; it is simply given as the 
average of tlie principal pressures pi {i = 1,2,3). 

The quantities Wi (i = 1,2,3), which arc the rcscaled principal pressures, encode the degree of 
anisotropy of the matter. U Wi = w Wi, the matter is isotropic. The density is a function of the 
metric and can thus be expressed as p = p{n, si, S2, S3). Equation (3) then entails 

dlogp 



w = 



dlogn 



Wi 



w ■ 



- 1 



9 log p 



dlogp 



(5a) 
(5b) 



^9 log Si ^3 dlog Sj 
We make the following simplifying assumption: 

Assumption 3. We suppose that the isotropic pressure and the density are proportional, i.e., we 
assume w = const, where 



we (-1,1). 



(6) 



According to Assumption 3, the density and the isotropic pressure behave like those of a perfect 
fluid with a linear equation of state satisfying the dominant energy condition. It seems natural 
to focus attention on anisotropic matter sources satisfying Assumption 3, since they generalize 
the behaviour of perfect fluid models widely used in cosmology. In the concluding remarks, see 
Section 8, we will indicate how to treat the more general case where w ^ const. 
Remark. In Assumption 3 the cases w = ±1 are excluded. The case w = 1 leads to different 
dynamics, which wc refrain from discussing here. The case w = —1 is excluded since it comprises 
as a subcase the dc Sitter spacctime which does not exhibit a singularity. 

Taking account of (5), Assumption 3 implies that 

p(n,si, 82,83) = n^+"''(/'(si, 52,53) , (7) 

for some function Tp{si, 82, S3). Moreover, for all i, the matter anisotropics Wi are functions of 
(si, 52,53) alone and 



■ Wi{Si,S2,S3) 



d log ip 
dlog Si 



9 log?/; 
d log Sj 



(8) 



(9) 



Finally we introduce a dimensionless time variable r by 

— = H-^ — - 

dr dt ' 

in the following a prime will denote differentiation w.r.t. t. 

Expressed in the new variables, the Einstein evolution equations decouple into the dimensional 
equations 

^ = -2x(l + ^^SfcEfc) , (10) 



H' = -3H 



1 --(!-«;) 



and an autonomous system of equations 



= -3n 

s' = -2si 



^{1 - w)T,i - {wi 



w) 



(i = 1,2,3) 
(i = 1,2,3) 



(Ua) 
(lib) 
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The Hamiltonian constraint results in 

O = 1 - where Y? := i (sf + + S^) , (12) 

which impUes < 1. This constraint is used to substitute for O in (11). 

The reduced dimensionless dynamical system (11) encodes the essential dynamics of Bianchi type I 
anisotropic spacetimes where the matter is described by the quantities Wi = Wi{si, 82,83) and 
w = const; once this system is solved, the decoupled equations (10) can be integrated and the 
solution (g'",A;'j) can be constructed. 

There exist useful auxiliary equations in connection with the system (11). In particular, the 
evolution equation for is given by 



ft' = n 



3{l-w)T,'' -Y^^Wk^k] . (13) 



Likewise, for p we get p' = —p [3(1 + w) + Wk S^]. 

The state space of the dimensionless dynamical system (11) is the four- dimensional bounded open 
connected set X given by 

X = {(El, S2, S3, 8US2, 8s) I (S^ < 1) A Sfc = 0) A (0 < Si < 1 Vi) A Sfc = l) } • 

This set can be written as the Cartesian product = 3C x T of two two-dimensional bounded open 
connected sets, 

3C={(Ei,S2,I]3)| (s2<l)A(^^Sfc=o)} . (14a) 

7={{su82,8s)\{0<8i<iyi)A(j2^Sk = l)} . (14b) 

The set % is the Kasner disc; it is typically depicted in a projection onto the plane with conormal 

(1, 1,1), see Figure 1. The boundary of 3C is the Kasner circle aOC = {(S^ = l)A(I]i+E2+E3 = 0)}. 
The Kasner circle contains six special points, which are referred to as LRS points: The three Taub 
points Ti, T2, T3 given by (Ei,E2,E3) = (2,-1,-1) and permutations, and the three non-flat 

LRS points Qi, Q2, Q3 given by (Ei, E2, E3) = (—2, 1, 1) and permutations. The six sectors of DC 
are denoted by permutations of the triple (123); by definition, E^ < Ej < E^ holds in sector (ijk). 

The set T is the interior of a triangle contained in the afiine plane si + S2 + S3 = 1, see Figure 1. 
The boundaries form a triangle with the three sides Si = 0, S2 = and S3 = 0; the corners are 
given by (1, 0, 0), (0, 1. 0) and (0, 0, 1). The functions 'ij){s\,82, S3) and Wi(si, S2, S3) in (7) and (8) 
are understood as functions on the domain T. 

We conclude this section with a discussion of the energy conditions. Our basic assumption is p > 0. 
Then the dominant energy condition is expressed in the new matter variables as \wi{si, S2, 53)! < 1, 
V(si,S2,S3) G T and Vi = 1,2,3; the weak energy condition is —1 < Vi = 1,2,3; the strong 
energy condition is satisfied if the weak energy condition holds and w > —1/3. We are mainly 
interested in matter models that satisfy the dominant energy condition; however we shall also 
discuss the more general case when the functions Wi are merely bounded: 

Assumption 4. We assume that the rescalcd principal pressures are bounded functions, 

sup [^^(si, S2, S3)| < 00 , Vi= 1,2,3. (15) 
T 



4 The flow on the boundaries of the state space 

The present section is devoted to a detailed analysis of the flow induced by (11) on the boundary 
of the state space. This analysis is rather technical but essential to understand the asymptotic 
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E2 S3 

(a) The Kasner disc % (b) The space T 

Figure 1: The four-dimensional state space X is the Cartesian product of the Kasner disc X and 
the set T. The latter is represented by (the interior of) a triangle; the center of the triangle is the 
point (si,S2,S3) = (1/3,1/3,1/3); for the sides, the values of (si, 82,53) are given in the figure, 
where the arrows denote the directions of increasing values (from to 1). 



behaviour (in particular toward the past) of solutions of (11). 

The state space X = XxT, see (14) and Figure 1, is relatively compact. Its boundary is the union 
of two three-dimensional compact sets, i.e., 

dX = (dXxT) U (^ X 07) . (16) 

The intersection of these boundary components is the two-dimensional compact set dX x 87. 

4.1 The vacuum boundairy: dOC X 7 

Consider the boundary component dX x T, which is characterized by = 1 (i.e., $1 = 0); we call 
this component the vacuum boundary — the reason for this terminology will become clear through 
the remarks in this subsection. The dynamical system (11) admits a regular extension from X to 
the vacuum boundary; we simply let — » 0, so that equation (11a) becomes = in this limit. 
This is independent of (si, S2, S3) € 7, since Wi is bounded on 7 by Assumption 4. The dynamical 
system (11) thus induces the system 



S' = , s' = -2si 



Si-V SfcSJ {i = 1,2,3), (17) 



on the vacuum boundary. 



The vacuum boundary dXx7is a solid torus whose cross section is 7; each cross section corresponds 
to (I]i,I]2,S3) = const and is an invariant subspace of (17). The fixed points of the system (17) 
are transversally hyperbolic and form a connected network of lines on dX x 87: 

• Kasner circles: There exist three circles of fixed points that can be interpreted as Kasner cir- 
cles. The Kasner circles are located at the vertices of T, i.e., let {ijk) be a cyclic permutation 
of (123), then KCj is given by 

KQ: S2 = l, {si,Sj,Sk) = {1,0,0). 
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• Taub lines: These are three hnes of fixed points given by 

TLi-. (Si, S,-, Sfc) = (2,-1,-1), {si,Sj,Sk) = {0,s,l-s),se[0,l]. 

• Non-flat LRS lines: These are three hnes of fixed points given by 

QL- (E„I]j,I]fc) = (-2,1,1), {si,Sj,Sk) = {0,s,l-s),se[0,l]. 

Let (ijk) be a cyclie permutation of (123). On the Kasner circle KCi there exist tliree 'Taub points': 
Tii, Tij, and T^^, which are given by (E„5]j-,Efc) = (2,-1,-1), (-1,2,-1), and (-1,-1,2), 
respectively. The point T^ (T^/j) is the point of intersection of KCj with TLj (TL^); Tu is an 
isolated Taub point since it does not lie on any of the Taub lines. Analogously, there exist three 
'non-flat LRS points' on KCi: Qu, Qy, and Q^fe, which are given by (Si,Ej, E^) — (—2,1,1), 
(1,-2,1), and (1,1,-2), respectively; the analogous comments apply. In Figure 2(a) we give a 
schematic depiction of the network of fixed points. 

Remark. Each of the fixed points on the Kasner circles is associated with a Kasner solution (Bianchi 
type 1 vacuum solution) of the Einstein equations, 

ds^ = -dt^ + oi t^P^dx^ ® dx^ + a2 f^^dx^ ® dx"^ + as t^^^dx^ ® dx^ , (18) 

where oi, 02, as are positive constants, and pi, p2, Pa are the so-called Kasner exponents, which 
are constants that satisfy pi + P2 + Pa = + P2 + P3 = 1- The relation between (Ei,E2,E3) 
and {pi,P2,P3) for a Kasner fixed point is 3pi = Ej + 1, i = 1,2,3. The fixed points on TLj 
and QLj represent the fiat LRS Kasner solution (Taub solution) and the non-fiat LRS Kasner 
solution, respectively; in these cases the Kasner exponents are (1,0,0) and permutations and 
(2/3,-1/3,-1/3) and permutations, respectively . 

Lemma 2. The solutions of the dynamical system (17) on the vacuum boundary are heteroclinic 
orbits, i.e., the a- and the ui-limit sets consist of one fixed point each; see Figure 2(b). 

Proof. Since Sj = const Vi, the result follows by studying the sign of s^ for each sector {ijk) (and 
at the special points) of the Kasner circle separately. □ 

Remark. Not only the Kasner fixed points themselves, but each solution of the dynamical sys- 
tem (17) on the vacuum boundary can be interpreted as a Kasner solution. This is simply because 
(Ei,E2,E3) = const and E^ = 1. 

4.2 The cylindrical boundary: 3C X 87 

The second component of dX is the set % x 9T; the set is a triangle and consists of the three 
sides {si = 0} {i = 1, 2, 3), hence 

5T=eiUe2Ue3, where := DC x {sj = 0} . (19) 

Each set Ci is compact; since it has the form of a cylinder, we call % x d7 the cylindrical boundary. 
By Ci we denote the interior of Cj. 

Asymptotic properties of matter models 

The dynamical system (11) admits a regular extension from X to the lateral surfaces of the cylinders 
Cj, cf. the discussion of the vacuum boundary; however, in general, the system is not extendible 
to the interior Cj, nor to the top/base surfaces. The simple reason for this is that, although 
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(a) Vacuum fixed point set 



(b) Dynamics on the vacuum boundary 



Figure 2: A schematic depiction of (a) the network of fixed points and (b) the flow of the dynamical 

system on the vacuum boundary dOC x T. Orbits in the interior of d% x T connect a transversally 
hyperbolic source with a transversally hyperbolic sink; for instance, the a-limit (w-limit) of the 
orbit through the point P e dX x T in Figure 2(b) is a point in sector (312) of KC2 (KC3). 



Wi{si, 82,83), i = 1,2,3, are bounded on T, cf. (15), in general these functions need not possess 
limits as (si, 32,53) converges to a point on 37. However, for reasonable matter models (such as 
collisionless matter and elastic matter), the matter anisotropics u;i(si, S2, S3) are functions on T 
that admit a unique extension to T. We make this a general assumption. 

Assumption 5. The matter anisotropies Wi{8\,82, S3) {i = 1, 2, 3) admit unique extensions from T 

to T. The extended functions are assumed to be sufficiently smooth on 87. (For the majority of 
our future purposes, continuity is a sufficient requirement; however, we will assume the functions 
to be to facilitate our analysis.) 

Assumption 5 ensures that the dynamical system (11) can be extended to C, (i = 1, 2, 3) and thus 
to the entire boundary dX of the state space X. The next assumption is a simplifying assumption 

on the functions Wi [i = 1, 2, 3) on 37. 

Assumption 6. We assume that Wi(si,S2,S3) = V- — constant on Ci (i.e., when Si = 0) for all 
i = 1,2,3. 

Rem,ark. Assumption 6 stems from basic physical considerations; in particular, it is satisfied for 
collisionless matter and for elastic matter, see Section 7. Loosely speaking. Assumption 6 means 
that the (rescaled) principal pressure in a direction i becomes independent of the values of gjj and 
gkk {i ^ j ^ k ^ i) in the limit ga 00 (<4> g" 0), provided that gjj, gkk remain bounded. (The 
complementary statement is implicit in our previous assumptions: The fact that Wi is well-defined 
for Sj = 1 means that the rescaled principal pressure in a direction i converges to a limit as gu — > 0, 
when gjj and gkk remain bounded from below.) 

For some matter models like elastic matter, there exists a function v{zi,Z2, Z3), which is defined on 
7 9 {zi, Z2, Z3), sufficiently smooth on 87, and symmetric in the arguments Z2 and 23, such that 

Wl(Sl, S2, S3) = 11(51, S2, S3) , «;2(si, S2. S3) = W(S2, S3, Si) , W3(S1, S2, S3) = W(S3, Si, S2) (20a) 

for all (si, S2, S3) € T. Equation (20a) follows from basic physical considerations; it is a consequence 
of the fact that there is no distinguished direction and reflects the freedom of permuting the axes. 
Assumption 6 is equivalent to the requirement that t;(0, Ci, C2) = f^- for all Ci, (2- 
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Definition 1. For matter models safAsfying the symmetry property (20a), in slight abuse of nota- 
tion we define a function v{s) on the interval [0, 1] by 



[0, 1] 9 s v{s) := v{s, 1 - 



s,0). 



(20b) 



The value of v at the endpoint s = is v{0) = V-, cf. Assumption 6; the value at s = 1 we denote 
by v+, i.e., v{l) = v+. 

For more general matter models, the initial data for the matter fields might break the symme- 
try (20a), an example being collisionless matter. Let us denote by / the initial data of the matter 
fields and by J the space of possible initial data. Let u be a permutation of the triple (123), which 
induces the permutation (si,S2,S3) '-^ (sct(i) j S(t(2) j So-(3)) on any ordered triple (si, 52,33) in M^; 
in slight abuse of notation we use the symbol a for the induced permutation as well. Finally, let 
/(ct) denote the initial data arising from / by the (induced) permutation. (For instance, in the case 
of collisionless matter, / represents the distribution function /o, and = /o o cr, cf. Section 7). 
Then there exists u:JxT— / x (^i, Z2, zs) ■u[/](2;i, Z2, zs), sufficiently smooth, such that 



where Ui is a permutation with <Ti{l) = i. Equation (20a') reduces to (20a) if there is no dependence 
on / S X; then w is simply replaced by v. (The symmetry of v in the second and third argument 
corresponds to the freedom of choosing even or odd permutations in (20a').) Assumption 6 is 
equivalent to the requirement that w[/](0, Ci, C2) = t^- for all I G I, for all (1, (2. For matter 
models of this type Definition 1 is replaced by 

Definition 1'. Fori = 1,2,3, let j, k be such that sgn{ijk) = +1. We define a function Ui{s) on 
the interval [0, 1] by 



The value of Ui at the endpoint s = is Ui{0) = V-, cf. Assumption 6; the value at s = 1 we denote 
byv+, i.e., Ui{l) =v+. Equation (20b') reduces to (20b) if there is no dependence on I and Ui is 
replaced by v. 

Remark. The important fact that Ui(l) is independent of z is a direct consequence of (20a') and 
the identity wi + W2 + W3 = 3w; see (21) and the discussion below. 

The dynamical system on Ci 

In the following let (ijk) be a cyclic permutation of (123). Wc begin our analysis of the flow of the 
dynamical system on the cylindrical boundary (19) by fixing the notation. Consider the cylinder 
Cj. Since Sj = on Cj, we may set sj = s, Sk = I — s, s G [0, 1]; the cylinder Cj is then given as 
the Cartesian product 3C x {0 < s < 1}. We call the set X x {s — 0} the base, and X x {s = 1} 
the top of the cylinder. The boundary of the top is the Kasner circle KCj , the boundary of the 
base is KC/j, see Figure 3. 

The matter anisotropies on the cylinder 6^ are represented by Wi, wj, and Wk- Assumption 6 implies 
that Wi = f_ on C^. When regarded as a function of s G [0,1], Wj coincides with the function 
Ui{s) of Definition 1', i.e., Wj{si,S2, S3) = '"[/(^-i)](sj, Sfe, Sj) = m[/(^-1)](s, 1 - s,0) = Ui{s) on Cj. 

Finally, Wk = 3w — Wi — wj, since Wi + Wj +Wk= 3w. Evaluation of this identity at s = 1 (i.e., 

for {si,Sj,Sk) = (0,1,0)) yields w_ + Ui{l) + = 3w, because Wk = V- for {si,Sj,Sk) = (0,1,0) 
by Assumption 6. Since, by definition, Ui{l) = V-^-, we obtain 



(20a') 




(20b') 



(21) 
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Using this equation in the form Sw — V- = V- + v+, the anisotropies on the cyUnder Cj take the 
form 

Wi=v- , Wj = Ui{s) , Wk = v++v-- Ui{s) , (22') 
where Ui{s) satisfies Ui(0) = V- and «,(!) = v+ (independently of i). 
In the simpler case (20a) we have 

Wi = V- , Wj = v{s) , Wk = v+ + V- — v{s) , (22) 

where v{s) satisfies v{0) — w_ and v{l) = v+, and again (21). However, by (20a), on C^, Wk can 
also be written as = v{l — s, 0, s) = v{l — s). Therefore, v{s) + v{l — s) = v^+v^; this identity 
can be written in an alternative form, 

Ka)-^)=-(Kl-)-^), (23) 

which states that v{s) — (f_ + f+)/2 is antisymmetric around s = 1/2. 

Using (22) we can write the dynamical system that is induced on by (11) as 

s' = -2s(l-s)(S, -Sfc), ^ = -Sn (^^{1 - w)Ei - [v- - w]^ , (24a) 

^ = -3fi Q(l - w)Ej - [v{s) - w]j , = -3n (^1(1 - w)Ek - [v+ +V-- v{s) - w]J , 

(24b) 

where wc recall that f2 = 1 — S^. The remainder of this section is devoted to a detailed analysis 
of the fiow of the dynamical system (24). 

In the case (22') the function v{s) in (24) is simply replaced by Uj(s). The qualitative dynamics 
of the resulting system remains unchanged because the key quantities in our analysis will turn out 
to be V- and v+, which are the same for (22) and (22'). 

In the following we analyze in detail the dynamical system (24) which is connected with the matter 
models satisfying (20a), (20b) and (22). The more general case, represented by (20a'), (20b') 
and (22'), leads to the same conclusions and will only be commented on sporadically. 

Flow on the lateral boundary of Gi 

The fiow on the lateral surface of Cj (which is given by = 1) is independent of the matter 
quantities, since = 0. Orbits on the lateral surface satisfy 

s' ^ <^ Sj- - Efe ^ , 

while s' = on the lines of fixed points TLj and QL^, see Figure 3. 

Flow on the bcise/top of Gi and definition of the various cases 

On the base of the cylinder Cj, where s = 0, we have Wi = wj = V- and Wk — v+, see (22) or (22'). 
The dynamical system induced on the base of 6, is thus 



i(l-w;)Sft-(w+-w) 



; (25) 



since 2v- + v+ = 3w we have v+ — w = 2{w — V-). 
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KG, 




Figure 3: The boundary component % x 87 consists of the three cylinders Ci U 62 U 63. In this 
figure, Ci is depicted together with the fiow of the dynamical system on the lateral boundary and 
the vacuum fixed points. 



Definition 2. We define the quantity f3 as 



/3:=2^ = ^. (26) 
1 — w 1 — w 



This quantity will play an essential role in our analysis. 

It is not difiicult to show that the solutions of the system (25) form a family of straight lines that 
are attracted by a common focal point, the fixed point 

Rfc: (Ei,I],-,Sfe) = /?(-!, -1,2), (27) 

see Figure 4. Analogously, on the top of the cylinder (which is given by Si = 0, Sj = 1, Sfc = 0) 
there is a fixed point Rj, where (Sj,5]j,Sfe) = /?(— 1,2,— 1). The flow of the dynamical system 
that is induced on the top of 6, is represented by a family of straight lines focusing at Rj. 

Depending on the value of f3, see (26), we distinguish several scenarios whose main characteristic 
is the position of the fixed point R^ in relation to the Kasner circle KCfe (on the base of C^). In 
the following we discuss these scenarios in detail, in particular in view of their compatibility with 
the energy conditions. 



C+ /3 > 1 (and /? < 2, see below). The point Kk lies beyond the Taub point on the Kasner circle 
(which is the point T^k on KCfe), i.e., S/sIr^ > 2; see Figure 4(a). This is the case if and 
only if v+ > 1, which is not compatible with the dominant energy condition. 

(] = 1. The point Hk coincides with the Taub point T^^, i.e., Sfc|Rs. = 2; see Figure 4(b). 
This is the case if and only if v+ = 1; accordingly, V- = {3w — l)/2 by (21), i.e., we have the 
chain {3w — l)/2 = V- < w < Vj^ = 1. Imposing the dominant energy condition on V- leads 
to the restriction 

At the same time the strong energy condition is satisfied. 
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A+ P G (0, 1). The point R^; lies between the center of the Kasner disk and the Taub point with 
a value SfelR^ G (0, 2); see Figure 4(c). This is the case if and only if {3w —1)/2<V- < w < 
v+ < 1. If w > —1/3, the dominant (and strong) energy condition are automatically satisfied 

for v± and /3 can take any value in (0, 1). If u> < —1/3, the condition u_ > —1 is stronger than 
the condition v- > {3w — l)/2; the range of (3 is then restricted to /3 G (0, 2(1 + ui)/(l — w)]. 

Ao /? = 0. The point lies at the center of the Kasner disc. This is the case iff v- = v+ = w. 

A_ (3 e (—1,0). The point Rfe lies between the center of the Kasner disk and the non-flat LRS 
point (the point Qkk) with a value Sfcl^j. G (—2,0); see Figure 4(d). This is the case iff 
— 1 + 2w < < w < V- < {1 + Vj)/2. If w > 0, the dominant (and strong) energy condition 
are automatically satisfied for v± and /? can take any value in (—1, 0). If u> < 0, the condition 
> — 1 is stronger than the condition v+ > — 1 + 2w] the range of (3 is then restricted to 

pe[-{i + w)/{i-w),Q). 

B_ (3 = —1. The point Rfe coincides with the non-flat LRS point Q^fe, i.e., ^kWk = ~2; see 
Figure 4(e). This is the case iff — 1 + 2w = v+ < w < V- = {1 + w) /2. The dominant energy 
condition requires 

u; > 

and at the same time the strong energy condition is satisfled. 

C_ (3 < —1 (and (3 > —2, see below). The point Rfe lies beyond the non-flat LRS point Qfefe, i.e., 
SfclRfc < —2; see Figure 4(f). This is the case if and only \i v+ < —1 + 2w < w < {1 + w) /2 < 
V-. To ensure compatibility with the dominant energy condition we must presuppose 

w>Q . 

The quantity /3 cannot assume arbitrary values less than —1 without a violation of the energy 
conditions. Imposing the energy conditions on v± restricts the possible range; we obtain that 
/? must be greater than or equal to the maximum max{— (1 + w)/(l — w), — 2}. The extremal 
case would be /? = —2, which corresponds to w_ = 1 v+ = 3w — 2). In this extremal case, 
the dominant energy condition requires w > 1/3. 

Remark. In a nutshell: Iff /? G [—2, 1], then compatibility with the energy conditions is possible; 
in general it is required that w be sufficiently large (w > 1/3 is always sufficient). 

Remark. For the function v{s) (where we recall v{0) = V- and v{l) ~ w+) we have the fundamental 
inequalities 

V- < w < v+ in the + cases {(3 > 0) 

V- = w = v+ in the case Aq {(3 = 0) (28) 
V- > w > v+ in the — cases (/? < 0) 

which will be used several times in the following. 

Remark. As we will see below, several characteristic properties of the flow of the dynamical system 
on the cylinder Cj change when |/?| > 2. Therefore, in the following, C+ will denote the case 
(3 G (1,2) only, and C_ will denote the case f3 G (—2,-1). The cases > 2 will be denoted by 
extra symbols: 

D+ (3 >2. The point Rfe lies far beyond the Taub point on the Kasner circle, i.e., SfelRj, > 4. 
Like case C+ this case is not compatible with the dominant energy condition. 

D_ (3 < —2. The point Rfe lies far beyond the non-flat LRS point Qfefe, i.e., Sfeln^ < —4. The only 
subcase that is compatible with the dominant energy condition is (3 = —2, which corresponds 
to V- = 1 (<^ v+ = 3w — 2). In this case, the dominant energy condition requires w > 1/3. 
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'J ^3 ^3 




(d) A_ ' (e) B_ ' (f) C_ 



Figure 4: Flow on the base of the cyhnder 6^; the orbits are straight hnes focusing at Rfc. The 
Kasner circle KCfc consists of fixed points that can act as sources, saddles, or sinks for orbits in 
the interior of the cylinder C^; in the figures, bold continuous [dashed] lines denote transversally 
hyperbolic sources [sinks]. In the cases B_|_ and B_, is not transversally hyperbolic. In the 
cases C_|_ and C_ the point is close to the Kasner circle since < 2; increasing the value 
of amounts to increasing the distance between R/j and the Kasner circle; in the case |/3| > 2 
the distance is greater than the Kasner radius. Accordingly, if /3 > 2 there appear transversally 
hyperbolic sinks in the lower half of the Kasner circle; if /3 < —2, some of the transversally 
hyperbolic sources in the upper half vanish. The flow on the top of the cylinder and the stability 
properties of the Kasner fixed points on KCj are obtained by a reflection w.r.t. the E^-axis. 



Stability of the fixed points on the boundary of Cj 

Combining the analysis of the flow on the lateral surface and on the base/top surface we obtain 
the stability properties of the flxed points on the Kasner circles KCj and KCfe. (We simply put 
Figures 3 and 4 on top of each other.) Depending on the case there exist (transversally hyperbolic) 
sources and sinks on KCj/KCt-; the location of these sources and sinks on KCfe is depicted in 
Figure 4 (the location of sources/sinks on KCj is obtained by a reflection w.r.t. the E^-axis). 
In cases B_|_ and B_, Rfe (which coincides with Tfefe and Qfcfc, respectively) is not transversally 
hyperbolic; the center manifold reduction theorem implies that, in case B_|_, the flxed point is a 
center saddle; in case B_, it acts as a sink, i.e., it is the w-limit for a two-parameter set of orbits 
in Ci ; the analog holds for Rj . 

The point Rj [Rfe] lies in the interior of the top [base] of the cylinder, iff we are in one of the A 
scenarios. By (24) we obtain 

(l-s)-i(l-5)'|j^^ =6/3, s-ifi'lj^^ =6/3. (29) 

In combination with the results on the flow on the base/top, see Figure 4, we obtain that Rj [Rfe] 
is a saddle in case A_(_, while it is a sink in case A_. 
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Remark. For completeness we briefly discuss the case Aq as well. In this case, /? = and the 
fixed points Kj/Kk are not hyperbolic. Center manifold analysis reveals that the character of 
these points depends on the derivative of v{s) at s = and s = 1, i.e., on v'{0) and v'{l). As 
the antisymmetry property (23) of the function v{s) guarantees that v'{0) = v'{l), there exist two 
subcases: 

Ao+ v'{0) = v'{l) < 0. Rj/R/s are center saddles. This case resembles the case A+. 
Aq— v'{0) = v'{l) > 0. Rj/Rfe act as sinks. This case resembles the case A_. 

Apart from a special subcase {v{s) = w), the degenerate case ^^'(O) = v'{l) = is excluded by 
Assumption 7 below. (If we analyze the system (24) with the function Ui{s) instead of v{s), 
cf. (22'), there might exist more subcases of Aq, because the signs of u^(0) and w^(l) might be 
different; however, since the focus of our analysis lies on the + and — cases and not on Aq , we 
refrain from discussing these additional subcases further.) 

The remaining fixed points on the boundary of Cj are located on the lateral boundary: TLj and 
QL^. Using (13), (21), and (26) we find 

n-^n' I^L. = 3(1 - V-) = 1(1 - w){2 + (3) , (30a) 
n-^n' Iqj^^ = 2(1 - v+) + (1 - V-) = 1(1 -w) + 1(1 - v+) = 1(1 - w)i2 - /?) . (30b) 

Therefore, if /3 < 2 each fixed point on QL^ acts as the source for exactly one orbit in the interior of 
Cj; if /? > 2, each point on QL- is the w-limit set for one interior orbit. The analog is true for TLj: 
If /3 > —2, then each fixed point on TL^ acts as the source for one interior orbit; if /? < —2, each 
fixed point on TLj attracts one interior orbit as r ^ oo. (The proof of these statements is based 
on the center manifold reduction theorem, where we use that the lateral boundary is the center 
manifold of QLj.) The borderline cases /? = ±2 cannot be dealt with by using local methods; these 
cases are discussed in Lemma 3 and Lemma 4 below. 



Dynamics in the interior of Gi 

The dynamical system (24) on the cylinder Cj admits a non-negative function, 

M(j) := (S j + /?) ' , M(',) = -3fi(l - «;)M(j) , (31) 
which is strictly monotonically decreasing on Cj whenever S, ^ —(3. The plane Sj = — /3 itself, i.e., 

= {(s,Si,S„Sfe)eej|l], = -/?} (32) 

is an invariant subset. It can be characterized as the unique plane, orthogonal to the base/top of 
the cylinder, whose boundary contains the points R^- and Kk- The plane CDj intersects the cylinder 
Cj whenever — 2 < /3 < 2. (When (3 = ±2, Dj is tangent to Cj; the intersection of the two is TLj 
or QLj. When |/3| > 2, Dj does not intersect Cj.) 

Lemma 3. Let 7 be an orbit in the (interior of the) cylinder Gi such that 7 ^ D,. The a-limit 
set of 7 is 

• one of the transversally hyperbolic sources on the Kasner circles KCj/KCfc (where the points 
on ODi are excluded), or 

• a fixed point on QL^ (when (3 <2), or 
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• a fixed point on TLj (when /? > —2). 

Each transversally hyperbolic source on KC^/KCfe is the a-limit set for a one-parameter family 
of orbits; each fixed point on QL^ and TLj is the a-limit set for one orbit in Cj (under the given 
assumption on (3). 

Proof. Assume first that < 2. The orbit 7 is either contained in the invariant subset = {S^ > 
— /?} or in the invariant subset §~ = {S^ < — /?} of 6^. Since the function M(^^-j is monotonically 
decreasing on the invariant sets §^ , the monotonicity principle imphcs that the a-limit set of 7 is 
contained on the boundary of however, is excluded, since Mj-j) attains its minimum there. 
This leaves the lateral boundary of 6^ as the only possible superset of the a-limit set of 7. From 
the structure of the flow on the boundary of 6^, which is depicted in Figure 3, we conclude that it is 
only the fixed points that come into question as possible a-limit sets. Combining these observations 
with the local analysis of the fixed points, which we performed in the previous subsection, leads 
to the statement of the lemma. In the case /? > 2, the set S^^ is empty and 8^ corresponds to the 
interior of the cylinder. The monotonicity principle, when applied to the function M(i), implies 
that 0(7) is located on the lateral boundary of Ci, where QL^ is excluded. The case /3 < —2 is 
analogous and thus the claim of the lemma is proved. □ 

Lemma 4. Let 7 C Cj, 7 ^ Then the u-limit set 0.1(7) 0/7 depends on the actual case under 
consideration (as characterized by the quantity /3): 

• (3 > 2 {-^ D^.). a;(7) is one of the transversally hyperbolic sinks on KCj/KCj. or a point on 
QLi- 

. /3 = 2 H D+). c.(7)CQL,. 

• /3g[-1,2) A, B, C+). o;(7)CD,. 

• (3 € (—2,-1) (<-^ C_). w(7) is one of the transversally hyperbolic sinks on KC^/KCfe or it 
lies on Di. 

• /3 < — 2 01(7) is one of the transversally hyperbolic sinks on KCj/KCfe or it lies on 

Remark. In analogy to the statement of Lemma 3, each transversally hyperbolic sink on KCj/KCfe 

is the w-limit set for a one-parameter set of orbits; for /3 > 2, each fixed point on QLj attracts 
exactly one interior orbit, so that the set QL^, when regarded as a whole, attracts a one-parameter 
set of interior orbits; see (30); in contrast, for (3 = 2, QL^ attracts every orbit in Cj, i.e., a two- 
parameter set. 

Proof. The proof is analogous to the proof of Lemma 3. Consider first the case < 2. On S^, 
which is the closure of the invariant subset defined in the proof of Lemma 3, the function M(i) 

attains its maximum on TL^; on §j , the maximum is attained on QL^; therefore, these sets are 
excluded as possible w-limits sets. Using the analysis of the flow on the boundary of the cylinder, 
the statement of the lemma ensues. The cases |/3| > 2 are analogous. □ 

Dynemiics on the invariant plane Dj 

Lemma 4 states that the future dynamics of orbits in Cj is connected with the properties of the 
flow on Dj. (This is for < 2, i.e., for all cases A, B, C.) We thus proceed by investigating the 
flow of the dynamical system on this invariant plane. 
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Let Di denote the set of fixed points on (the interior of) for the dynamical system (24). The 
coordinates (s, Sj, S/j) of any fixed point P G Dj must satisfy T,j = Sfe = /?/2 and 

v{s) = ^iv.+v+), (33) 

where s G (0, 1). In the + and — cases (i.e., /3 ^ 0) we have V- ^ v^; therefore there exists at least 
one solution of equation (33). (Analogously, in the cases Ao+ and Aq—, existence of a solution is 
guaranteed by the fact that the derivatives of v{s), or Ui{s), at s = and s = 1 have the same 
sign.) 

Remark. Because v{s) — {v- + v+)/2 is antisymmetric around s = 1/2, cf. (23), one solution of 

equation (33) is s = 1/2. Consequently, the point P^^'' with coordinates (s, Sj, fife) = (1/2)(1, (3, (3) 
is always a fixed point on D^. (If we analyze the system (24) with the function Ui{s) instead of 
v{s), cf. (22'), s = 1/2 will in general not be a solution.) 

We make the simplifying assumption that the zeros of the function v{s) — (w_ + v+)/2 are simple 
zeros: 

Assumption 7. We assume that the function v{s) either satisfies 

v'{s) ^0, for all s e [0, 1] such that v{s) = ^~ , (34) 
or v{s) = w for all s. 

Assumption 7 is automatically satisfied for collisionless matter and clastic matter and helps to 
avoid unnecessary clutter in the analysis of the flow on D,; requiring (34) is equivalent to assuming 
that the set of fixed points on Di is a discrete set of hyperbolic fixed points. In addition, it 
follows that the number of fixed points is always odd, i.e., = 2d + 1, d G N. Accordingly, we 
can write = {P^ . . . ^\ pf'^ , P^^"* , . . . , P^*^' } , where the coordinates of these fixed points 
satisfy < < ... < < S(o) < < ... < S(d) < 1- The only exception we admit 

is the case v{s) = w, which is a subcase of Aq ; in this special case, the set Dj of fixed points is 
not a discrete set but coincides with the line Ej = = in Dj; this special case will be briefly 
commented on at the end of this section. 

Remark. The antisymmetry property (23) of v(s) guarantees that s„ + s_„ = 1, n = 0, 1, . . . , d. 
This is not true in general, if we analyze the system (24) with the function Ui{s) instead of v{s), 
cf. (22'). 

Let P e Dj be a fixed point in Dj with coordinates (s, /3/2, /3/2). The eigenvalues of the linearization 
of the dynamical system (24) at P are given by 



-A(i_^)(2+/3)(2-/3) 



1 , L 645(1 - s) 

- 3(1-^)^(2-/9) -^^-) 



(35a) 



and 

-^(l-«^)f^|p- (35b) 

The former are associated with eigenvectors tangential to the plane Dj, the latter is associated 
with an eigenvector transversal to Dj. It follows that P is a sink if v'{s) > 0, while P is a saddle 
if v'{s) < 0. Therefore, by Assumption 7, the set Dj = {P^"'^-', . . . ,Pj''^} of fixed points can be 
regarded as an alternating sequence of saddles and sinks. In the + cases, inequality (28) implies 
the structure Dj = {sink, saddle, sink, . . ., sink, saddle, sink}; in particular, d+1 fixed points are 
sinks and d points are saddles. In the — cases, inequality (28) implies the structure Dj = {saddle, 
sink, saddle, . . ., saddle, sink, saddle}; in particular, d+1 fixed points are saddles and d points are 
sinks. 
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Remark. In the cases A-|- the fixed points Rj and Rfc can be added to the alternating sequence 
{Pj- P^''''} of saddles and sinks; this follows from the local analysis of the fixed points Rj/R^. 

Accordingly, in the case A+, the sequence {Rfc, P-~'^\ . . . , P-''^ Rj} is of the type {saddle, sink, 
saddle, . . ., saddle, sink, saddle}, while in the case A_, it is of the type {sink, saddle, sink, . . ., sink, 
saddle, sink}. The cases Ao+ and Ao_ can be subsumed under the cases A+ and A_, respectively, 
the only difference being that Rj/R/c are non- hyperbolic. 

The next step in our study of the flow of the dynamical system on the plane Dj is to investigate 
the global dynamics. To this end we consider the function 



N = {1 - Y;'^)-\{s) , 
where k{s) is positive and satisfies the differential equation 



sil-s)- 



ds 



(36) 



(37) 



By (28), in the + cases the function k{s) goes to oo as s and s ^ 1; in the — cases, k{s) in 
the limit s — > and s — > 1. A straightforward computation shows that N is strictly monotonically 
decreasing on Dj except on the fixed point set Dj, i.e.. 



N' 



w 



N . 



and 



N" 



■Y^(2 + /3)'(2- 



{v{s) 



l-s)- 



(38) 



N . 



The monotone function N on excludes the existence of periodic orbits, homoclinic orbits and 
heteroclinic cycles in (the interior of) Dj. Using the monotone function in combination with the 
local analysis of the fixed points, it is possible to solve the global dynamics of the flow on Dj. Since 
the relevant arguments are rather standard (because the theory of planar flows is well-developed) 
we omit the derivation here and merely summarize the results in Figure 5. Note that in the cases 
B+, C+, the flow on dDi is given by a heteroclinic cycle; in the case B+ it is represented as follows: 



^33 



T 



kk 



Having solved the dynamics on Dj we are now in a position to complete the analysis of the future 
dynamics of orbits in Cj. 



Dynamics in the interior of Gi (cont.) 

While Lemma 3 describes the past asymptotics of orbits in Ci\Di completely, the description of 
the future asymptotics by Lemma 4 is incomplete in the cases where \(3\ < 2. Exploiting the results 
of the previous subsection we are able to complete the discussion. 

Lemma 5. Consider either of the + cases (with [3 < 2), i.e., either of the cases A_|_, B4-, C+ 
(or Ao+j. Let 7 be an orbit such that ^ (jL1)i. Then 0^(7) is one of the fixed points of Dj. More 
specifically, of the #Di = 2d+l fixed points, there are d+1 sinks that attract a two-parameter set 
of orbits each, while each of the d saddles attracts a one-parameter family. 



19 



(a) B+, C+ 



(b) A_,_ (and Ao+) 



(c) A_ (and Ao_) (d) B_, C_ 



Figure 5: The figures display the flow of the dynamical system on the plane D^. We assume < 2 
(cases A, B, C) so that Dj is an invariant subset lying in Cj. The number of fixed points in the 
interior of CDj is necessarily odd, i.e., #Dj = 2d+l, d G N; we depict the case d= 1, i.e., #Di = 3. 
[The case =f^Di = 1, which is not depicted, is commented on in square brackets.] Non-gcnoric orbits 
are represented by dashed lines, (a) The boundary dT>i forms a hcteroclinic cycle in cases B_|_ and 
C-)-. Generic orbits converge to this hcteroclinic cycle as r ^ — oo. The central fixed point P^"'' is 
a saddle; the two other points P^^^' are sinks. [In the case #Dj = 1, is the a-limit set, and 

the (unique) fixed point vf*^ is the cj-limit set for all (non-trivial) solutions.] (b) In case A-|_, the 
properties of the three fixed points Dj are analogous. The points Rj/Rfe are saddles; one interior 
orbit converges to Rj, one to Rfe as r ^ — oo. (Whether these two orbits converge to the sinks 
as r ^ oo, as depicted, or the to central fixed point, depends on the function v{s).) The past 
attractor for generic orbits consists of the Kasner points in the upper left/lower right corner. [In 
the case #Di = 1, the fixed point pI""* is the w-limit set for all (non-trivial) solutions.] The fiow in 
the case Ao+ resembles the fiow in the case A+; however, Rj/Rfc are center saddles, (c) In case 
A_, the central fixed point P^""* is a sink, the two other interior fixed points P^^"'^'' arc saddles; 
the points Kj/Kk are sinks. The Kasner points in the upper left/lower right corner are the a- limit 
for generic orbits. [In the case #Dj = 1, the fixed point Pf^^ is a saddle, which leaves Rj/Rk as 
the exclusive future attractor for generic orbits.] The flow in the case Aq- resembles the flow in 
the case A_; however, Rj/Rk are not hyperbolic, (d) The cases B_ and C_ are similar to A_; 
however, the Kasner points in the lower left/upper right corner are sinks. [In the case #Di — 1, 
the fixed point Pf^^ is a saddle, which leaves these Kasner points as the generic future attractor.] 



Proof. The flow on is known in detail from the analysis of the preceding subsection. We see that 
the only cases that require a careful analysis are B+ and C-|-, since in these cases, the boundary 
dTfi of the plane T>i is a hcteroclinic cycle, which has to be excluded as a possible oj-limit of 7. To 
this end we use the function N given by (36). Suppose that 7 converges to the hcteroclinic cycle 
dDi. A straightforward computation shows that 



1 



w 



-(I 



«;)(E|-/?2) 



(39) 



along 7. The term in square brackets on the r.h.s. of (39) is a positive function, which does not 

converge to zero along 7; on the contrary, for sufficiently large times t this function is approximately 
equal to (3/2)(4 — /3^) for most of the time (which is because the orbit 7 spends most of its time 
in a neighborhood of the fixed points on dDi). In contrast, the second term on the r.h.s. of (39) 
converges to zero as r ^ 00. It follows that TV ^ along 7. However, this is a contradiction, 
because N is infinite on dDf, therefore, the hcteroclinic cycle dD^ is excluded as a w-limit of 7. □ 
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(d) A_ (and Ao_) (e) B_ (f) C_ 



Figure 6: A schematic depiction of the flow on the cyhnder Cj. In contrast to Figure 5, to simphfy 
the illustration we assume in this figure that the function v{s) is such that there is only one fixed 
point Dj in the plane Di. 



Lemma 6. Consider either of the cases A_, B_ (or Aq- j. Let 7 be an orbit such that 7 ^ Di. 

Then uj^j) is one of the fixed points on Di, i.e., one of the points Rj, Rfc, or one of the points 
of the set . More specifically, of the #Di = 2d+ I fixed points, there are d sinks that attract a 
two-parameter set of orbits each, while each of the d+1 saddles attracts a one-parameter family. 
The fixed points Rj, Rk (which coincide with Qjj, Qkk in B-) attract a two-parameter family of 
orbits each. 

Lemma 7. Consider the case C_. Let 7 be an orbit such thai 7 ct 2),;. Then uj{^/) is one of 
the fixed points of the set or one of the transversally hyperbolic sinks on the Kasner circles 
KCj/KCfc. More specifically, of the f/=Di = 2d + I fixed points, there are d sinks that attract a 
two-parameter set of orbits each, while each of the d+1 saddles attracts a one-parameter family. 
Each sink on the Kasner circles attracts a one-parameter family of orbits. 

Proof. The statements of the two previous lemmas follow from Lemma 4, the analysis of the flow 
on Di, and the local analysis of the fixed points. □ 

Lemma 3 describes the past dynamics, Lemma 4 (in the cases D) and Lemmas 5-7 (in the cases A, 
B, C) describe the future dynamics of orbits in Cj. The qualitative behaviour of orbits is depicted 
in Figure 6. 

Remark. In Assumption 7 we allow for a special case, defined by the requirement v{s) = w for all 
s. In this special case, is not a discrete set of fixed points but a line of fixed points, the line 
= Sfc. It follows from (35) that the fixed points of Dj are transversally hyperbolic sinks. Using 
the monotone function (36) and the simple structure of the flow on Dj it is not difficult to show 
that a;(7) e Dj for aU orbits 7 C ej\Dj. 
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5 Local dynamics 



In this section we perform a stability analysis of the fixed points in the (closure of the) state space 
X. This local dynamical systems analysis is an essential ingredient in order to understand the 
global dynamics of solutions in X. 



Fixed points on the boundary dX 

In the previous section we have analyzed in detail the dynamical system (11) on dX. The fixed 
points of this system reside on the cylindrical boundary Ci U 62 U 63 C dX. Since Ci is given as 
Sj = 0, the direction orthogonal to Cj is given by the variable Sj. Prom (11) we obtain 

L,=o = -2 i^i - -'^^ - (1 - . (40) 

where we again use the convention {sj,Sk) = (s, 1 — s). Evaluated at the fixed points we get 

QLi: sr's'i\^^^^ = 6, (41a) 

(41b) 



TLi : 


s,- ^s'A „ = -6, 











KC.VKCfc : sri^: ^^^^ = -2(S, - S^y^) , (41c) 

where (3 is given by (26); (3 encodes the dependence on the rescaled matter anisotropics, which are 
represented by v±, cf. (22). The results of (41) can be combined with the rcsidts of Section 4 to 
complete the stability analysis of the fixed points on the boundary of the state space. The results 
are summarized in Table 1, Table 2, and Table 3. 

Example. Take, for instance, the fixed point in the case A+: Prom Pigure 4(c) we see that Rfc 
is a sink on the base of the cylinder C^; sec also Figure 6(c). Figure 5(b) shows that there exists 
exactly one orbit in Dj (and thus in Ci) that emerges from R^. Together with (41c) it thus follows 
that Rfc possesses a two-dimensional unstable manifold (which lies in the interior of the state space 
X), so that there exists a one-parameter family of orbits (in X) that converge to Rfc as r ^ —00; 
cf. Table 1. 



Fixed points in the interior of X 

The dynamical system (11) possesses a set of fixed points in the interior oi X (= DC x T), which is 
given as the set of solutions of the equations 

Si = S2 = S3 = , w\ = W2 = w^i = w; (42) 

the first condition defines a point in %\ in general, the second condition defines a subset of T. Since 
the principal pressures coincide, each of these fixed points represents the fiat isotropic Friedmann- 
Robertson- Walker (PRW) perfect fluid solution associated with the equation of state p = wp. 

Assumption 8. Por simplicity, we assume that the matter model is such that w\ = = = w 
defines a point in T; hence (42) defines an isolated fixed point in X , which we denote by F. 

Remark. This is an assumption out of convenience. (It is automatically satisfied for the matter 
models we consider in detail in Section 7 — coUisionless matter and elastic matter.) The analysis 
of this section is equally valid when (42) has more than one solution. In that case, P denotes any 
of the multiple fixed points in X . 
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point eCiC dX 



point G QLj 
d+1 points e Di 
d points G Di 
Rj/Rfc 
point P G KCj/KCfc 
[condition on P] 
point G TLi 



B. 



1-parameter family 

one orbit 
1-parameter family 

1- parameter family 

2- parameter family 

[1 < Sj/fc|p] 



1-parameter family 
one orbit 

1- parameter family 

part of dDi 

2- parameter family 
[1 < E,./fc|p < 2] 



1-parameter family 
one orbit 

1- parameter family 

2- parameter family 
[1 < S,./fc|p < 2/P] 



does not attract any interior orbit as r ^ — oo 



Table 1: In the + cases, the fixed points on dX are a-limits for orbits of the interior of the 
state space X; in this table we list how many interior orbits converge to a particular fixed point 
as r — oo. The fixed points on the Kasner circles that satisfy the condition in the table are 
transversally hyperbolic sources and thus act as the a-limit for a two-parameter set of orbits each. 
(The points that do not fulfill the condition do not attract any orbits.) The fixed points on QL^ 
(where the end points and Q^^ are excluded) act as the a-limit for a one-parameter set of 
interior orbits; hence, in total, a two-parameter set of orbits has an a-limit on QL^. In case B+, 
the fixed points Rj/Rfe (which coincide with Tjj/Tkk) he on 92) j, which is a heteroclinic cycle. 
The collection of these cycles, dDi U U dD^, forms a so-called heteroclinic network. 



point € Gi C dX 




A_ 


B_ 


C_ 


point G QL^ 


a 


1-parameter family 


1-parameter family 


1-parameter family 


point P G KCj /KCfe 


a 


2-parameter family 


2-parameter family 


2-parameter family 


[condition on P] 




[1 < S,/fe|p] 


[KS.aIp] 


[l<S,/fc|p] 


point G TLj 




does not attract any interior orbit as r 


—^ —00 or r — * 00 


d+1 points G Dj 


U) 


2-parameter family 


2-parameter family 


2-parameter family 


d points G Di 


U) 


3-paramctcr family 


3-paranic;ter family 


3-parameter family 


Rj/Rk 


U) 


3-parameter family 


3-parameter family 




point P G KCj /KCfe 


U) 


does not act 


as w-limit 


2-parameter family 


[condition on P] 




[- 


■] 


[S,./fe|p < 2/(3] 



Table 2: In the — cases, some fixed points on dX are a-limits, some are w-limits for orbits of the 

interior of the state space X; in this table wc list how many interior orbits converge to a particular 
fixed point as r ^ —00 (denoted by a) or r ^ 00 (denoted by w). In case B_, the fixed points 
Rj/Rk coincide with Qjj/Qkk- The proof that these points attract a three-parameter set of orbits 
requires center manifold analysis. 
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point eGiC dX 




D+ 


D_ 


point e QLj 




not attractive 


1-paramcter family 


point P e KCj/KCfe 


a 


not attractive 


2-parameter family 


[condition on P] 




H 


[^<^j/k\p] 


point e TLj 


a; 


not attractive 


1-paramcter family 


point P e KCj/KCfe 




not attractive 


2-parameter family 


[condition on P] 




[-] 


[5^j/fe p < -1] 



Table 3: If > 2, the properties of the fixed points change considerably. If /3 > 2 (D+), neither 
of the fixed points can be the a/w-limit set of an interior orbit. (However, a fixed point can be 

contained in the a-limit set of some orbit 7 C X.) If /3 < — 2 (D_), then each fixed point P 
on KCj/KCfe with 1 < Sj//c|p is a transversally hyperbolic source, while Sj//c|p < —1 yields a 
transversally hyperbolic sink. 



To study the stability properties of F, it is useful to define a different set of coordinates on 1. Let 
{ijk) be a cyclic permutation of the triple (123) and set 



1 + + e*. 



Si = 



1 + e** + e*J' 



Sk = 



(43) 



This defines a coordinate system 9 (ti,tj) T, such that the point (ti,tj) = (0,0) corresponds 



to the center (si, S2, S3) = (1/3, 1/3, 1/3) of J; the limit Sj ^ corresponds to U 



-00, Si 



Oto 



tj —00; the limit s, ^ 1 corresponds to U ^ 00, Sj ^ 1 to — > 00; finally, s/s — > corresponds 
to a combined limit ti ^ 00, tj ^ 00. In these coordinates the fixed point F is given by 



F : (Si,S2,S3)|F = (0,0,0), {ti,tj)iF := 



Using {ti,tj) on T, the rescaled pressures (8) can be expressed as 



Wj = w + 2 — log^p , 



Wk 



2lt 



dti dt 



(44) 



(45) 



where ip is regarded as a function of the two variables (t^, tj). Since w\ = = W3 at F, it follows 

that (ti, fj) G T is a critical point oi tp; wc set tp = tj). 

The dynamical system (11) can be represented in the four variables (S^, T,j, ti, tj) as 
"1 



Si 



-30 



-(1 - u;)S„ - {Wn - w) 



4 = -2(2S„ + S„), (n,m) e {{i,j),{j,i)} 



where Sl = l— |(Sf + T,'j + SjEj). The matrix representing the linearization of (11) at the point 
F contains the Hessian f) of V' evaluated at (ti,ij). The four eigenvalues of this matrix are given 
by 
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A±± = 7(l-H -1± Vl-A± 



where A± are given by A± = 64^ 



tr f) + f)i,- ± V(tr f) + f)i,)2 - 3 det f) /[3(l-w;)2]. 



Lemma 8. // (i) \) is positive definite, i.e., detl) > and trf} > 0, then all eigenvalues X±± 
have negative real part and F is a hyperbolic sink. If (ii) () is negative definite, i.e., dett) > and 
tr [) < 0, then the eigenvalues A+± are positive and A_± negative so that F is a hyperbolic saddle 
with a two- dimensional stable and a two-dim,ensional unstable manifold. If (iii) dot f) < 0, then 

A±+ and A have negative real part, whereas \+- has positive real part, so that F is a hyperbolic 

saddle with a three-dimensional stable and one- dimensional unstable manifold. In the exceptional 
case (iv) det f) = 0, there exists at least one zero eigenvalue. 
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Proof. From the inequality 

detf,<^trf)2-f)|., (4g) 

we obtain (tr f) + f)jj)^ — 3det f) > 0, hence A± are real numbers. Studying the sign of A± in the 

different cases, the claim of the lemma follows immediately. □ 



6 Global dynamics 



On the state space X we define the positive function 

M=(1-E2)-i^(si,s2,S3), (47) 

where ip{si, S2, ss), (si, 32,53) G T, is defined by (7). By a straightforward computation, where 
we use that M can be written as M = 3H^/n^^^ (which makes it possible to apply (10) and the 
equation n' = —3n, which follows from d\/det g/dt = 3H\/det g), we obtain 

M' = -3(1-w)S2M. (48a) 

The computation of higher derivatives reveals that 



M" 



S2=0 



= -9(1 - w)M ^{wk - wf (48b) 



on the subset = of the state space. Therefore, since w < 1 by Assumption 3, we obtain 
M' < when 7^ and M"'|s2=o < except at the point F, cf. Assumption 8. Hence we 
have proved that M is a strictly monotonically decreasing function along the flow of the dynamical 
system (11) on A'\F. 

Lemma 9. Let 7 be an orbit in X\F. Then the a-limit set 0(7) and the iv-limit set 0^(7) satisfy 

• either 0(7) = F or a{j) C dX, and 

• either u{'y) = F or uj{j) C Ci U 62 U 63 C dX, cf (19). 



Proof. The lemma is a direct consequence of the monotonicity principle. Since M is strictly 
monotone along the flow of the dynamical system (11) on A'\F, the a/w-limit set of 7 must 
be contained in F U dX. Because these limit sets are necessarily connected, we find that either 
a(7) = F or a(7) C dX, and likewise for 0^(7). Moreover, M\gxxT = 00 (i.e., M ^ 00 along every 
sequence converging to a point on dOC x T), hence this subset of dX is excluded for 0^(7), which 
leaves uj{-f) C 5C x = Ci U 62 U 63. □ 



It is the quantity /3, see (26), that determines the details of the asymptotics. Therefore we discuss 
the cases /3 > and (3 < separately; the former includes the + cases A+, B+, C+, and the 
extreme case D+, the latter includes the — cases A_, B_, C_, and the extreme case D_. 



The + cases 

Theorem 1 (Future asymptotics). In the + cases A+, B+, C+, and D+, the fixed point F is 
the future attractor of the dynamical system (11); every orbit in the state space X converges to the 
fixed point F as r — > +00. 
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Interpretation of Theorem 1. Since the fixed point F corresponds to a FRW perfect fluid solution 
associated with the equation of state p = wp, the theorem states that each Bianchi type I model 
with anisotropic matter that satisfies V- < w < v+, cf. (28), isotropizes toward the future and 
behaves, to first order, like an (infinitely diluted) isotropic perfect fluid solution. 

Remark. Theorem 1 can be proved under assumptions that are considerably weaker than the ones 
we made. (Assumptions 5-7 are rather irrelevant for the future asymptotics — they are tailored 
to the past asymptotics.) This fact will become obvious from the proof and is thus not further 
commented on. 



Proof. To prove Theorem 1 we must show that ^(7) C Ci U 62 U C3 is impossible, see Lemma 9. To 

this end we use that V'(si, S2, ^3) ^ 00 as (,si, S2, .S3) 87 (which wc prove below). It then follows 
that M = sup M = 00 on the cylindrical boundary and consequently on the entire boundary dX 
of the state space, cf. (47). To finish the proof we apply the monotonicity principle [22] for the 
function M on X\{F}: The monotonicity principle precludes the possibility that 6^(7) might be 
contained on dX for any orbit 7 C A" and the theorem is established. To show that ~* 00 as 
87 is approached, we use that V- <w <v+, which holds by assumption, cf. (28). Employing the 
coordinates (43) on T and the representation (45) of the rescaled matter quantities, we find that 
there exists some e > such that for all tj there exists T > and 

8 Wi-w 8 Wi-w 

-logV = ^-<-e, _log^ = ^->e 

for all ti < — \T\ and ti > T respectively. Consequently, \ogip 00 and thus ^ ^ 00 as ti ^ ±00. 
(Likewise, ip ^ 00 as tj ^ ±00.) We thus conclude that i/' = 00 on 87 and the theorem is 
proved. □ 

Theorem 2 (Past asymptotics). Lei 7 C X\F be an orbit of the dynamical system (11). 

• In the + cases A+, B+, C+, the a-limit set 0/7 is 

A+ one of the fixed points on 8X; 

B+ one of the fixed points on 8X or, possibly, the heteroclinic network 8D1 U dT>2 U 8T>3; 
C+ one of the fixed points on 8X or one of the heteroclinic cycles 8T)i, 8D2, 8Dz. 

Convergence to a fixed point is the generic scenario, i.e., the set of orbits that converge to 
the heteroclinic structures is a set of measure zero in X. Table 1 gives a complete list of the 
fixed points on 8X together with the number of orbits that converge to the particular points. 

• In the case D+, the a-limit set of 7 is represented either by a heteroclinic cycle or by a 
heteroclinic sequence on 8X and thus contains a sequence of Kasner points; in particular, 
a(7) cannot be a fixed point. 



Interpretation of Theorem 2 and conclusions. Theorem 2 (in conjunction with the results of Table 1) 
entails that in each of the + cases, generic solutions approach the fixed points on the Kasner 
circle(s), which implies that the past attractor of the dynamical system (11) is located on the 
Kasner circle(s). Therefore, each generic Bianchi type I model with anisotropic matter, where the 
matter is supposed to satisfy the requirements of the + cases, approaches a Kasner solution toward 
the singularity. However, there are interesting differences among the cases. In case A^. each Kasner 
solution is a possible past asymptotic state, i.e., given an arbitrary Kasner solution (18), represented 
by the Kasner exponents {pi,P2,Pi) = (l/3)(I]i + l, S2 + 1,^3 + 1), Pi +P2+P3 = P1+P2+P3 = 1) 
there exist anisotropic Bianchi type I solutions that converge to this Kasner solution as t ^ 0. 
In particular, each of the three Taub solutions, which are characterized by {pi,pj,pk) = (1,0,0), 
represents a possible past asymptotic state. In case B+ this is no longer correct. Since the Taub 



26 



points are (center) saddles there does not exist any solution that converges to a Taub solution 
as t ^ 0. However, the Taub points are part of a larger structure, the heteroclinic network 
dDi U U 9D3, which represents a potential a-limit set for (a probably two-parameter set) 
of solutions. Accordingly, we conjecture that there exists a two-parameter set of solutions whose 
asymptotic behaviour is characterized by oscillations between different representations of the Taub 
solution. In case C+ (which is characterized by 1 < /? < 2) there is an entire one-parameter set 
of Kasner solutions (including the Taub states) that are excluded as a-limit sets. The condition 
of Table 1 yields that a Kasner solution is excluded as a past asymptotic state for Bianchi type I 
models with anisotropic matter iff max; > 2//3, or, equivalently, max-ipi > (2 -|- /3)/(3/3). In 
other words, generic Bianchi type I models behave asymptotically like Kasner solutions as the 
singularity is approached; however, the set of possible Kasner limits is restricted to those that are 
"sufficiently different" from the Taub solutions. If /3 > 2, i.e., in case D+, this set is the empty set. 
Each point on the Kasner circles takes the role of a saddle point. This induces generic oscillatory 
asymptotic behaviour of Bianchi type I models. The a-limit of a generic solution contains a 
sequence of Kasner states, hence the solution undergoes a sequence of phases ("Kasner epochs"), 
in each of which the behaviour of the solution is approximately described by a Kasner solution, 
and "transitions" between these epochs. 

Remark. Theorem 2 (in conjunction with the results of Table 1) also implies that there exist 
(non-generic) Bianchi type I models whose asymptotic behaviour is not connected to the dynamics 
of any Kasner solution. There exist solutions whose asymptotic behaviour is characterized by 
(Si, Ej, Efc) /3(— 1, 1/2, 1/2) as the singularity is approached, and solutions with (Ej, Ej, Efe) 
/3(2,— 1,— 1), where the latter concerns the case A+ since /? < 1 is required. Furthermore, in 
the cases B+ and C-|-, there exist solutions whose asymptotic behaviour is oscillatory: An orbit 
approaching the heteroclinic cycles/network {dDi, dD2, dD^) corresponds to a solution that is 
represented by a periodic sequence of "Kasner epochs" associated with the Kasner fixed points that 
lie on the cycles; in the Kasner epochs the solution is characterized by (Ej,Ej,Efc) = (— /3, /3/2± 
V374^/4^, (3/2 T x/a/iyi^) . 

Remark. The distinction into the different + scenarios is intimately connected with the energy 

conditions; see the list of cases in Section 4. If the energy conditions hold (case A_|_, /3 < 1), each 
Kasner solution represents a possible past asymptotic state for anisotropic Bianchi type I models. 
If the energy conditions are satisfied only marginally (case B+, (3 = 1), i.e., at the onset of energy 
condition violation, the Taub solutions are excluded as past asymptotic states. Finally, in the case 
/3 > 1 (case C-j- and case D+) the energy conditions are violated; in this context, the quantity /3 — 1 
can be regarded as a measure for the magnitude of the violation. Increasing the value of /3 — 1 from 
zero to one leads to an ever increasing exclusion of Kasner states from the past attractor; if /3 — 1 
reaches one (/? = 2), the set of Kasner states that represent past asymptotic states has shrunk to 
the empty set. This is the onset of oscillatory behaviour of generic solutions. 



Proof. In the proof we choose to restrict ourselves to the cases A+, B+, C+, since these are the 
difficult cases; the proof in the case D+ is similar (but simpler). The first step to prove Theorem 2 

is to show that a (7) must be contained in dX. Using Lemma 9 this amounts to proving that 
a(7) = F is impossible. Therefore, assume that there exists 7 C X\F such that a(7) = F. 
By Theorem 1 we know that 0^(7) = F, hence 7 must be a homoclinic orbit. However, this a 
contradiction to the fact that the function M, which is well-defined in X and in particular at F, 
cf. (47), is strictly monotonically decreasing along 7. Consequently, a(7) = F is impossible, and 
we conclude a(7) C dX. 

In the next step of the proof we show that 0(7) is contained in a subset of the boundary dX: 
0(7) C (aCi U 962 U dQa) U (Di U D2 U D3) . By definition, dX is represented by the disjoint union 

dx = (aeiuae2u5e3)u(DiuD2UD3)u((eiue2ue3)\(DiuD2UD3)) u{dXx7) , 

ldX]i [ax]2 ldX]3 [9 AT] 4 
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sec Section 4. Suppose that there exists an orbit 7 such that a{'-f) 9 P G [SA'Ja; w.l.o.g. P G 6^ 
(P ^ Dj) for some i. Since the a-hmit set a{'j) is invariant under the flow of the dynamical system, 
the entire orbit through P, jp, and its w-limit must be contained in 0(7) as well. By Lemma 5, 
tj(7p) is a point Pd^ of the set D;. However, Pd^ is a saddle in X, see Table 1; therefore, since 
Pd; G 0^(7)) either 7 is contained in the unstable manifold of Pd^, or the intersection of 7 with 
the unstable manifold is empty. The first alternative contradicts the assumption P e 0(7). In the 
second case wc appeal to the Hartman-Grobman theorem to infer that 0(7) contains not only Pd^ 
but also points of the unstable manifold of Pd^- (Consider a sequence t„ such that t„ —00 and 
l{Tn) Pd; as n ^ 00. By the Hartman-Grobman theorem there exists a sequence f„ = Tn + Srn, 
with Sr-n > 0, f„ < T„_i, such that the sequence 7(t„) possesses an accumulation point on the 
unstable manifold of Pd;-) The local dynamical systems analysis, cf. Table 1, thus implies that 
a(7) contains points that do not lie on Ci but in the interior of X. This, however, is a contradiction 
to Lemma 9. The assumption has led to a contradiction; therefore, the intersection of a(7) and 
[9 A"] 3 must be empty. The procedure to show that a (7) fl [9 A"] 4 = is identical, where we use the 
flow on dX X T, see Figure 2(b), and the stability analysis of the fixed points as summarized in 
Table 1. We have therefore shown that 0(7) C (SCi U 8^2 U ^63) U (Di U D2 U D3) . 

In the last step of the proof we investigate which structures on the set dCi U 9C2 U 9C3 are 
potential a-limit sets. Since this set is two-dimensional, possible a-limit sets are fixed points, 
periodic orbits, homoclinic orbits, or heteroclinic cycles/networks. The analysis of Section 4, see 
Figures 3, 4, and Figure 6, shows that, in case A_|_, there do not exist any other potential a-limit 
sets on dCi U 9C2 U 9C3 than fixed points. (This concludes the proof of the theorem for case A+.) 
In case B+, there exists an entangled network of heteroclinic cycles, a 'heteroclinic network'; in 
case C_|_, there are three independent heteroclinic cycles: dDi, dT>2, dD^ To complete the proof 
of the theorem we have to investigate whether a(7) can coincide with dT>i for some i in case C-|- 
(where 1< /3 < 2). 

Consider the heteroclinic cycle dDi for some i. At the point P on dTii given by (sj, Sj,Sk) = (0, 0, 1), 
Si + /3 = 0, Tij = Efe, we consider the three-dimensional hyperplane Tip : Sj = Sfe orthogonal 
to dT>i; the natural axes on Tip are Si, Sj, and Si -I- f3. We study the Poincare map induced on 
Hp by the dynamical system. The subspace Sj = is an invariant subspace of Hp (since 6^ is an 
invariant subspace of X). Equation (31) and Figure 5(a) reveal that P is a saddle point for the 
dynamical system on Hp, the unstable manifold being the s^-axis, the stable manifold being the 
(Si + /3)-axis. 

A straightforward calculation shows that 

srh', = 30-V3^/I^X, (Ae[-l,l]) 

along dDi. Under the assumptions of case C+ the r. h. side is strictly positive for all A, hence, for 
some e > 0, s~^Si > 2e along dDi and s^^s'i > e in a sufficiently small neighborhood Ui of dDi 
in X. Therefore, as long as a solution is contained in this neighborhood Ui of dDi, we find that 
Si{T) is bounded by some constant times e^^ (for decreasing t). We infer that a solution either 
remains within Hi for all sufficiently small r (t ^ —00) or |Si + /3| becomes large (and the solution 
leaves Ui in this way). Letting d<; = — (3/2)(l — 'w)^dT (so that decreasing r-time corresponds to 
increasing ^-time) we find that 

^ (Si + /3) = (Si + /3) - - V-) , (49) 

where Wi = Wi{si,S2, S3) satisfies — > V- as Sj 0. 

Let us consider a simplified problem. A differential equation of the type /'(c) = /(?) + g{<i) with 
a function (7(c) that converges to zero as c ^ 00 possesses a unique solution such that /(<;) ^ as 
C — > 00; if, initially, / is larger than this special solution, then /(<;)—> 00 as ^ —> 00; if, initially, / 
is smaller, then /(c) — > —00 as c — > 00. 
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Likewise, we obtain from (49) that the quantity S; + /? is increasing if it is positive and sufficiently 
large initially, while it is decreasing if it is negative and sufficiently small initially. In the former 
cases, the solution eventually leaves Uj, since Sj + /? has become too large; in the latter case, the 
solution leaves Hi , because T,i + (3 has become too small. Consider a one-parameter set of initial 
data, e.g., Sj = const, sj = const' in Tip. Let (Sj + /3)o be an initial value of Sj + /? such that the 
associated solution satisfies (Sj +/?)(?) > 5 for sufficiently large ? (for some appropriately chosen 6). 
By continuous dependence on initial data, there exists an open interval of initial data containing 
(Sj + P)o such that the analog holds for each initial datum of this interval. Let (Sj + /3)pos denote 
the infimum of (Ej + /3) of the maximally extended open interval. (This infimum exists because 
there exists an analogous maximally extended open interval comprising the small initial values of 
Sj-|-/3 which lead to + < —6.) By construction, the solution with initial datum (Sj + /3)pos 
must remain in Ui for all times, which implies that the a-limit of this solution is dDi. Accordingly, 
dDi is a possible a-limit set of orbits in X 

To sec that convergence to 91); is non-generic (in the sense stated in the theorem) we can invoke 
the Hartman-Grobman theorem for discrete flows (for the flow on Hp); we obtain that dDi has 
the character of a saddle: There exists solutions that converge to dDi, but generically, solutions 
are driven away from this cycle. Alternatively, we note that in order to converge to dDi a solution 
must satisfy 

(Si + /3)(0) = (Si + /3)o = / e-' [wi{si{c;),Sj{c;),Sk{<;)) - v.] d<; , 

which is a direct consequence of (49). Since (Sj + /3)o = o(l) as Sj ^ (and uniformly in Sj), 
there cannot exist an open neighborhood of P in Hp such that solutions of the Poincare map with 
initial data of that neighborhood possess P (i.e., dDi) as an a-limit set. This completes the proof 
of the theorem. □ 



The — cases 

To simplify the analysis of the — cases, we require the matter model to satisfy a typical 'generic- 
ity' assumption which concerns the neighborhood of the isotropic state of the matter: We assume 
definiteness of the Hessian f) of V' at the Priedmann point (i.e., we exclude the cases (iii)-(iv) of 
Lemma 8). The matter models we discuss in Section 7 satisfy the required assumption automati- 
cally. 

Lemma 10. In the — cases A_, B_, C_, and D_ the fixed point F is a hyperbolic saddle with 
a two-dimensional stable and a two-dimensional unstable manifold. 

Proof. In terms of the coordinates (43) on T, the fixed point F is given by (Si, S2, S3) = (0, 0, 0) 
and {ti,tj) = {ii,ij), where {ii,tj) is the (single) critical point of the (positive) function ip on 
T, see (44). Using that tjj = on dl (which we prove below) it follows that (ti,tj) is a global 
maximum of tp. Therefore, at least for a matter model with a generic characteristic function tp, 
the Hessian i) of tp at {ii,ij) is negative definite. Part (ii) of Lemma 8 then implies the statement 
we have wished to prove. It remains to show that as 9T is approached. To this end we use 
the representation (45) of the rescaled matter quantities and the fact that V- > w > v+, which 
holds by assumption, cf. (28). Then there exists some e > such that for all tj there exists T > 
and 

d Wi-w d Wi-w 

-logtP = —->e, _logV>=^-<-e 

for all ti < — |r| and ti > T respectively. Consequently, log^ — > —00 and thus V' — * as t j — > ±00. 
Hence ip = on. dl. □ 
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Theorem 3 (Future and past asymptotics). Let 7 C X\F be an orbit of the dynamical 
system (11). In the — cases A_, B_, C_, and D_, one of the following possibilities occurs: 

• 0(7) is a fixed point on dX and uj{'^) is a fixed point on dX; this is the generic case; 

• a{'^) is a fixed point on dX and ^(7) — F; 

• a(7) = F and w(7) is a fixed point on dX. 

Table 2 gives a complete list of the fixed points on dX together with the number of orbits that 
converge to the particular points. 

Proof. The theorem is essentially proved by repeating the first two steps of the proof of Theorem 2. 
The main difference is the character of the fixed point F, which is described in Lemma 10. Note 
that on 9Ci U 8^2 U dQ^ (which is the relevant part of dX where a{'^)/io{"f) must reside) there 
do not exist any structures (like periodic orbits or heteroclinic cycles) that qualify as possible 
a-/a;-limits except the fixed points of Table 2. □ 

Interpretation of Theorem 3 and conclusions. Theorem 3 (in conjunction with the results of Table 2) 
entails that in each of the — cases, generic solutions converge to fixed points on the Kasner circle(s) 
as r — > — 00, which implies that the past attractor of the dynamical system (11) is located on the 
Kasner circle(s). Therefore, each generic Bianchi type I model with anisotropic matter, where the 
matter is supposed to satisfy the requirements of the — cases, approaches a Kasner solution toward 
the singularity. In contrast, toward the future, the asymptotic behaviour of generic solutions is 
quite diverse. In particular, the asymptotic properties of the matter model, in the shape of the 
function v{s), play an important role, because these determine the number of fixed points in 
Dj. Assuming the simplest case, #Dj = 1, which is in accord with the examples of anisotropic 
matter models discussed in Section 7, we obtain the following results: In the cases A_ and B_ 
the behaviour of generic solutions towards the future is (S^, E^, Tik) f3{2, —1, —1). In case C_, 
generic solutions converge to a Kasner solution towards the future. However, there is only a 
subset of Kasner solutions that qualify as possible future asymptotic states; it is only those Kasner 
solutions with max; S; > 2//3 (which corresponds to max; pi > {2 + /3)/(3/3) when written in terms 
of the Kasner exponents) that come into question. This set of Kasner solutions that qualify as 
future asymptotic states becomes larger with decreasing (3. Finally, in the case D_(/3 < —2), 
generic solutions converge to Kasner solutions toward the future, and conversely, each Kasner 
solution occurs as a future asymptotic state. 

Remark. Theorem 2 (in conjunction with the results of Table 2) also implies that there exist (non- 
generic) Bianchi type I models whose asymptotic behaviour is quite different. In particular, in 
cases A_, B_, C_, the orbits converging to as r ^ 00 give rise to solutions whose asymptotic 
behaviour toward the future is characterized by (S^, Sj , Sfe) /?(— 1, 1/2, 1/2). If #Dj = 1, this 
behaviour is non-generic; however, for matter models whose properties are such that #Dj > 1, 
this asymptotic behaviour is shared by a generic set of solutions. The most interesting non-generic 
behaviour concerns isotropization: There exist non-generic solutions (a one-parameter set) that 
isotropize toward the past, and non-generic solutions (a one-parameter set) that isotropize toward 
the future. 

Remark. In Section 4 we have seen that the — cases A_, B_, C_ (and the special case /3 = — 2 of 
case D_) are compatible with the energy conditions. However, despite this fact the assumptions 
of the — cases are representative of a matter model that is rather unconventional. This is due 
to the fact that the isotropic state of the matter is not energetically favorable and thus unstable 
(cf. the proof of Lemma 10, where we showed that the Hessian f) of V' is negative definite at the 
coordinates of F). 
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(a) C+ (1< /3 < 2) 



(b) B+ (/3 = 1) 



(c) A+ (0 < /3 < 1) 




(d) A_ (-K /3 < 0) (e) B_ (/3 = -1) (f) C_ (-2 < (3 < -1) 



Figure 7: A schematic depiction of the projection onto the Kasner disc of interior orbits. By 
abuse of notation, we use the same letter to denote a fixed point and its projection onto the 
Kasner disc. The center of the disc is the projection of the fixed point F (which represents the 
isotropic solution). Generic orbits are in bold. In the cases B_|_/C_|_, the dashed lines represent the 
heteroclinic cycles/network. In the case B_|_, the heteroclinic network could also attract orbits as 
T — > — oo (not depicted). In the case A_, there are Kasner states which act both as a- and w-limit 
of interior orbits. 



The case Aq 

The case Ao represents a bifurcation between the + cases and the — cases. The analysis of 
this case is difficult for several reasons, the most obvious being the fact that center manifold 
analysis becomes ubiquitous, since several fixed points possess large center subspaces. Another 
problem is the fact that the asymptotic behaviour of the function ip{si, S2, S3) as (si, S2, S3) 07 
is undetermined by the assumptions of case Aq ; there are several alternatives: ~^ 0, "0 — > cxd, -0 
might converge to some positive number, or -0 might not be convergent at all as (si, S2, S3) — *■ 87. 

Our expectation is that — under consistent assumptions on the behaviour of the function on T — 
the subcase Ao+ resembles the case A_|_, and the subcase Aq- resembles the case A_. However, 
a careful analysis is required to establish these vague expectations as facts; this goes beyond the 
scope of the present paper. 
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7 Anisotropic matter models 



In this section we present in detail three important examples of matter models to which one 
can apply the main results of this paper: Collisionless matter, described by the Vlasov equation, 
elastic matter, and magnetic fields. For a general introduction to collisionless matter and the 
Vlasov-Einstein system we refer to [1, 17]; the Bianchi type I case is discussed in detail in [6]. For 
a thorough discussion on the general relativistic theory of elasticity we refer to [2, 4, 9, 10, 20]; we 
will confine ourselves to deriving the energy-momentum tensor of elastic bodies. 



7.1 Collisionless matter 



Consider an ensemble of particles with mass m that move along the geodesies of a spacetimc (M, g) . 
(The geodesic motion of the particles reflects the condition of absence of any interactions other 
than gravity; in particular, collisions are excluded.) The ensemble of particles is represented by a 
distribution fimction ('phase space density') / > 0, which is defined on the mass shell, i.e., on the 
subset of the tangent bundle given by g{v,v) = —m?, where v denotes the (future directed) four 
momentum. If (t, a;') is a system of coordinates on M such that dt is timelike and d^i is spacelike, 
then the spatial coordinates of the four momentum are coordinates on the mass shell, and we 
can regard / as a function / = f{t,x'^,v^), i,j = 1,2,3. The distribution function / satisfies the 
Vlasov equation 

dtf+ ^d,j f - l^T^^^vVd^i / = , (50) 

where w° > is determined in terms of the metric g^i, and the spatial coordinates of the 
momentum via the mass shell relation g^vV^v'^ = —m?. The energy-momentum tensor is defined 
as 

(51) 



In the definition of a Bianchi type I solution of the Einstein- Vlasov system, where we use the 
notation and conventions of Section 2, / is assumed to be independent of a;' and so (50) takes the 
form 

dtf + 2k\v'd,,f = 0. (52) 
It is well known, see, e.g., [12, 13, 16, 14], that the general solution of (52) is given by 

f = f{t,v^) = foivj) , (53) 

where /o is an arbitrary function which corresponds to the initial data; usually, /o is assumed to 

be compactly supported, so that the integral (51) is well defined. Using the spatial metric gij, 
cf. (1), we may replace dot g by det g in (51) and we can write |uo| = y^rn^ + gijV^vi . 

Accordingly, equation (51) means that the energy- momentum tensor is given as a function of the 
spatial metric gij (which depends on the initial data /o); we have thus verified Assumption 1 for 
collisionless matter (in Bianchi type I). If /o satisfies the symmetry requirement 

fo{vi,V2,Vs) = fo{-Vl,-V2,V3) = fo{-Vl, V2, -V3) = fo{vi, -V2, -V3) , 

and the metric is diagonal, then the energy-momentum tensor T^^^ is diagonal as well, which proves 
that Assumption 2 is satisfied. Therefore, the Bianchi type I Einstein- Vlasov system admits the 
class of diagonal models as solutions, cf. the remarks at the end of Section 2. 

The rescaled principal pressures are given by 



Wi = 



fovl 



-1/2 



dvidv2dv3 



s,. I fa vf Efe Skvl 



-1/2 



dvidv2dv3 



1/2 



dv\dv2dv3 



1/2 



(54) 



dv\dv2dv3 
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when we assume that the particles have zero mass (see, however, the remark at the end of this 
section). As a consequence, Assumption 3 holds with 



w = 



(55) 



and (7) takes the form 

k 

Taking the limit Sj ^ in (54) we find that 

_ Sj J fov'^lsjv'^ + Skvl]-'^/^dvidv2dv3 



1/2 



dvidv2dv3 . 



(56) 



Wi = V- = . 



I fo[sjVj + Skvl]'^/'^dvidv2dv3 



Wk = 1- W-i 



(57) 



on Cj, i.e., Assumptions 4-6 are satisfied (with v- = 0). In order to make contact with (20a') we 
identiiy I with /o and define 



U[I]{ZI,Z2,Z3) : = 



Zl I foivi,'"2,V3)vf Y^k^kvl 



-1/2 



dvidv2dv3 



j fo{vuV2,Vz) Y^k^kvl 



1/2 



(58) 



dvidv2dvs 



for all /o and (2:1,22,-23) S T. Using that /(^r) = /o o cr it is straightforward to show that (20a') 
holds, i.e., Wi = o o-j leads to (54). The function Ui{s) of Definition 1' can be read off 

easily; equation (57) thus corresponds to (22') with v- = 0, v+ = 1; cf. also (21). We conclude 
that /3 = 1, i.e., an ensemble of massless coUisionless particles constitutes an anisotropic matter 
model that is of class B^.. 

Finally we note that Assumption 7is satisfied as well. One can show that Ui{s) (which replaces v{s) 
in (34)) is strictly monotonically increasing, hence there exists only one solution of (33) and thus 
only one fixed point on Di, i.e., #Dj = 1. Analogously, one can convince oneself that Assumption 8 
holds. 

Remark. If wc consider ensembles of coUisionless particles with positive mass, m > 0, Wi and w 
are functions of (si, S2, S3) and an additional scale, which can be taken to be x or n = {detg)~^^^ 
(where we recall that x = n^/^(siS2S3)~^^^)- (In [6], x is replaced hj z = m^/{m^ + x).) In 
Section 8 we will show that the analysis of this paper carries over straightforwardly to this more 
general situation (the main reason being that the length scale is a monotone function). 



7.2 Magnetic fields 

For an electromagnetic field represented by the antisymmetric electromagnetic field tensor F^v the 
energy-momentum tensor is given by 

T^M = _ J_ ( F/i Z?« _ ix;* pa ] 

47r V 4" 1/-^ a-^ /3y • 

The equations for the field are the Maxwell equations 

V^F'^" = , V^^F^^] = . 
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Consider specifically a purely magnetic field in a Bianchi type I spacctime with metric (1) that is 
aligned along, say, the third axis. In this case the electromagnetic field tensor takes the form 

/O 0\ 
_ K 
~ -if ' 

\o oy 

where K determines the magnetic field: = 0, = 0, B"^ = K{g^^g'^'^g^^Y/'^ . For a diagonal 
metric, the Maxwell equations imply that is a constant; hence the energy density p, 

p = 1- , 

is a function of the metric (which depends on the initial data for the magnetic field). Furthermore, 
is diagonal and 

T\ =T\=p, T\ = -p . (59) 

Accordingly, Assumptions 1 and 2 are satisfied and diagonal models exist. It follows from (59) 
that wi = 1, W2 = 1, ■if 3 = —1, and w = 1/3. For these normalized anisotropic pressures all the 
remaining assumptions are satisfied straightforwardly and we have V- = 1, v+ = 1 and thus 

/3 = -2. 

Therefore, in our classification, the asymptotic behaviour is of type D_. 

The conclusions are identical if we add an electric field parallel to the magnetic field. The Maxwell 
equations show that = E'^ = and E^ = L{g^^g'^'^g^^)^^'^, where L = const. The energy density 
becomes Snp ~ g^^g^^(K^ + L^), the energy-momentum tensor remains a functional of the metric, 

and (59) and its consequences remain valid. 

Remark. The magnetic field considered in [11] is not aligned with one axis and thus not included in 
our analysis. In fact, when not aligned, the magnetic field has to rotate and its dynamics becomes 
non-trivial and has to be added to the system of equations; in addition, the energy-momentum 
tensor is no longer diagonal. 

7.3 Elastic matter 

In elasticity theory, an elastic material in a completely relaxed state is represented by a three- 
dimensional Riemannian manifold (iV, 7), the material space, where the points of N identify the 
particles of the material (in the continuum limit) and 7 measures the distance between the particles 
in the completely relaxed state; let us denote by X^, A = 1, 2, 3, a system of local coordinates on 
N. The coordinates on the spacetime (M, g), on the other hand, are denoted hy x^, /j. = 0, . . . ,3. 
The state of the elastic material is described by the configuration function which is a (smooth) 
map 

such that the kernel of the deformation gradient Tip : TM — > TN is generated by a (future-directed 
unit) timelike vector field u, i.e., keiTtjj = (u) or u^d^ip^ = 0. The vector field u is the matter 
four-velocity; by construction, ■^~^(p) (i.e., the world-line of the particle p € N) is an integral 
curve. 

The pull-back of the material metric by the map ip is the relativistic strain tensor 

h^u = d^i>^d„ip^ JAB ; 
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since hfj^^u^ = 0, h represents a tensor in the orthogonal eomplement {u)^ of u in TM. h^^i, is a 
Riemannian metric on {u)^ , hence /i^j, has three positive eigenvalues hi, /12, h^- Note also that 
t^uh^v = 0, i.e., h is constant along the matter flow. 

The material is unstrained at the point x iff hf^^v = g^v holds x, where 

is the Riemannian metric induced by the spacetime metric g on (w)"*". The scalar quantity 

n = ^det gh = \/hih2hz (60) 

is the particle density of the material. This interpretation is justified by virtue of the continuity 
equation 

{nu^) = . 

Definition 3. An elastic material in a Bianchi type I spacetime is said to be Bianchi type I sym- 
metric iff C^-hf^i, = 0, where ^i, i = 1,2,3, are the Killing vectors of the spacetime. Furthermore, 
if the matter four velocity u is orthogonal to the surfaces of homogeneity, the elastic material is 
said to be non-tilted. 



According to Definition 3 (and Cuh^v = 0), the strain tensor of a Bianchi type I symmetric and 
non-tilted elastic material satisfies 

/ioo = hok = , hij = const . (61) 

in the coordinates of (1). Consequently, /i'^- = g^'^hjk and thus hi, /12, /13 depend only on t. 

A specific choice of elastic material is made by postulating a constitutive equation, i.e., the func- 
tional dependence of the (rest frame) energy density p of the material on the configuration map, 
the deformation gradient and the spacetime metric. An important class of materials is the one for 
which this functional dependence enters only through the principal invariants of the strain tensor. 
In this case we have 

P = p(qi,q2,q3) , (62) 

where 

qi=tr/i, q2=tr/i^, q3 = tr/i^; 

since = (qf — 3qiq2+2q3)/6, one of the principal invariants can be replaced by the particle density 
n. The materials described by (62) generalize the class of isotropic, homogeneous, hyperelastic 
materials from the classical theory of elasticity. In the following we shall refer to these material 
simply as elastic materials. The stress-energy tensor associated with these materials is obtained 
as the variation with respect to the spacetime metric of the matter action Sm = — J ^/\g\ P- The 
general expression which results for the stress-energy tensor is 

T,. = 2^-pg,., (63) 

which results in (3) for the metric (1). 

In Bianchi type I symmetry, and assuming that the elastic material is non-tilted, it follows from (61) 
that the principal invariants of the strain are functions of the spatial metric g only, hence non-tilted 
elastic materials satisfy Assumption 1. In the following we introduce a natural class of constitutive 
equations that are compatible with Assumptions 3-8 as well. 

Remark. There exist elastic materials which do not satisfy these assumptions. For instance, for the 
constitutive equation considered in [3], which was taken from [10], the rescaled matter quantities 
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(written in terms of the variables si, 52,53) do not have a continuous extension to the boundary 
87 of the space T. Therefore, strictly speaking, the analysis of this paper does not apply to 
this class of materials. However, the reason for this problem is that the variables (si, 52,53) are 
ill-adapted to this particular constitutive equation; replacing the 'triangle variables' (,si, S2, 53) by 
'hexagon variables', cf. [3], remedies this deficiency and Assumptions 3-8 hold w.r.t. this alternative 
formulation of the problem. The dynamics of the elastic matter models investigated in [3] seems 
considerably more complicated than the dynamics of the matter models considered in this paper; 
however, modulo the necessary reformulation of the problem, it appears that the models of [3] are 
of type 13 = 2. 

We consider elastic materials whose constitutive equation depends only on the number density n 
and the dimensionless shear scalar s defined by 



„2/3 



Whuh2,h3) J 



(64) 



where A and S denote the arithmetic and geometric mean functions, respectively. The fundamental 
inequality .71(01,02,03) > 3(01,02,03), which holds for all real non-negative numbers 01,02,03, 
with equality iff oi = 02 = 03, see [5], implies that s > and s = iff /ii = /i2 = h^. Thus the 
shear scalar s measures deviations from isotropy. Using the requirement p — p{n, s) in (63) leads 
to 

f dp \ 2 dp f, 1 ^ 

T^, = pUf,u^ + l^n— -pj g^, + \ h^.u - ■^^iQ^.u 

We restrict ourselves to constitutive equations of the quasi Hookean form 

p = p{n) + jl{n) f{s) , (65) 

where p{n) is the unsheared energy density and /i(n) the modulus of rigidity. For such a constitutive 
equation we obtain 

Tij.t. = pu^u„+pg^u + /'(s) ^/i^^ - ^qig^u^ . (66) 

Here, p denotes the isotropic pressure p, which is given by 

p = p{n) + i>{n)f{s), p = n^^(^), ,y = n^^(^). (67) 



dn \ n I ' dn 



The quantity p{n) is the unsheared pressure. Note that if the constitutive equation (65) depends 
only on n, i.e., if /(s) = 0, then the elastic material reduces to a perfect fluid with energy density 
p and pressure p. 

To specify the functions p and /i in the constitutive equation (65) we postulate a linear equation 
of state between the unsheared pressure and the unsheared energy density, i.e., p = ap , and a 
linear equation of state between the modulus of rigidity and the unsheared pressure, i.e., jj. = bp. 
By (67) this is equivalent to setting 

p = pon"'^^ , jl = poab n"+^ 

for some constant po > 0. Accordingly, the density p and the pressure p take the form 

p = pon"+^(l + o6/(s)) , and p = ap (68) 

and the stress-energy tensor (66) becomes 



T 



r^/n \ ah f'is) ( hnv 1 



(69) 



We make a number of assumptions on the parameters a, b and the function /: 
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(Al) a e [-1,1), ab> 0. 

(A2) /(s) > and /(s) = if and only if s = 0. 
(A3) f'is) > 0. 

Conditions (Al) and (A2) imply that the energy density is non- negative and has a global minimum 
at zero shear. Condition (A3) then corresponds to the (physically reasonable) assumption that the 
rest energy of the body is an increasing function of the shear. 

Remark. If 6 = 0, the modulus of rigidity /t vanishes and the elastic matter reduces to a perfect 
fluid with linear equation of state p = ap. If a = (so that p = 0), the choice of b is irrelevant, 
since ab = 0; this is clear because shear cannot occur for dust. 



Denoting by Wi the rescaled principal pressures, equation (69) leads to 

«,, = ^=a + 2(3 + s)-^^(^-i) . (70) 
p 'l + abf{s) Vqi 3/ 

^ V ' 

=: Q{s) 

(A4) We assume that sup^^g Q{^) < ^^'^ that the limit Qoo •= l™2^cxi Q{z) exists. 

For clastic matter that is Bianchi type I symmetric and non-tilted the relativistic strain is given 
by (61). We assume that /i^ (which plays the role of the initial data for the matter) is diagonal; 
by scaling the spatial coordinates we can achieve hij = Sij. Substituting into (69) it follows 
immediately that T*^(to) is diagonal (since gij{t()) is diagonal), hence T^j remains diagonal for all 
times by the evolution equations (2a). Therefore, the class of diagonal Bianchi type I models with 
elastic matter is well-deflned, cf. the remarks at the end of Section 2. 

Using hij — Sij we further obtain 

hi=g^\ /l2=^/'^ /^3=/^ qi=x, n= {deig)-^'^ , (71a) 

and therefore 



s 



Inserting (71) into (70) we find 



3 = (siS2S3)"'/^-3. (71b) 



Wi{suS2,S3)=a + 2(^Si-^Q{s) (72) 

and w = a. This shows that the elastic materials under consideration satisfy Assumption 3. 
Furthermore, by (A4), the functions Wi can be extended to 87; on Cj we obtain 

Wi = a- ^(3cx> , = a + 2 - Qoo , Wk = a + 2 ^Sfc - ^ Qoo ; 

hence Assumptions 4-6 hold as well, where w_ = a — 2Q^/3. Evidently, the rescaled pressures are 
of the form (20a) with 

^^(^1, Z2, Z3) = a + 2 ^zi - Qoo ; 
the function v{s) is given by v{s) = a + 2{s — l/3)Qao, see (20b); cf. also (22). We find that 

v-=a - |Qoo , v+ = a+ ^Qoo , and /? = ^^^^^ ■ (73) 
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(A1)-(A4) imply that /? > 0, i.e., clastic models arc anisotropic models of the + type (A_|_, B_|_, 
C+, or D+) or of the Aq type. However, the latter case occurs only if Q^o = 0, which implies 
that v{s) = w (= a), which is the special subcase of Aq introduced in Assumption 7. It is easy to 
sec from the linearity of the function v{s) that Assumption 7 is satisfied for all values of Qoo- In 
particular, there exists only one solution of (33) and thus only one fixed point on B^, i.e., = 1. 

Remark. Assumption (A3) excludes the possibility Qoo < 0, hence neither of the — cases can occur 
for the elastic materials. A formal way to obtain an elastic material that is of type — would be to 
allow /' < 0; this is unphysical, however, since it implies that the rest energy of the elastic body 
has a maximum at zero shear (instead of a minimum). 

It remains to check the validity of Assumption 8. By (72), the equations Wi = a are satisfied 
only if si = S2 = S3 = 1/3 (i.e., at the center of T) or at points (si, 52,^3) where Q{s) = 0. 
By (A3), the latter equation has no solutions for finite values of s (i.e., in T), hence Assumption 8 
holds. 

Remark. It is immediate from (72) that the dominant energy condition is satisfied if 
(A5) sup^>o <3('2) < -4(a) :=niin{|(l-a),|(l + a)}. 

Clearly, < A(a) < 1; A{a) = if and only if a = ±1, while A{a) = 1 if and only if a = -1/3. If 
in addition a, > —1/3 holds, then the strong energy condition is also satisfied; note in particular 
that the validity of the strong energy condition does not impose any condition on the constant b. 

Example (John materials), /(s) = s. This class of elastic materials was introduced in [20]. In this 
case Qoo = 1) hence /3 = 4/[3(l — a)] by (73). Consequently, this class of materials is of type A+ if 
a e [-1,-1/3), of type B+ if a = -1/3, of type C+ if a e (-1/3,1/3), and of type D+ if a > 1/3. 
Condition (A5) is satisfied if and only if 

a = -i , -1 < 6 < , (74) 

i.e., in case B_|_. 

Example (Logarithmic materials), /(s) = log(£ + s) — loge with e > 0; for simphcity, £ = 3. Since 
Qoo = 0, these materials are of type Aq where v{s) = w {= a). Condition (A5) reads 

ab < A{a) , (75) 

hence these materials satisfy the dominant energy condition if b is sufficiently small (a ±1)- 

Example (Power-Law materials), /(s) = (e + s)^ — with e > and A > 0; for simplicity, 
£ = 3. In this case Qoo = \ hence ^ = 4A/[3(1 — a)] by (73). Consequently, this class of 
materials is of one of the types A_|_, B_|_, C-|_, or D_(_, depending on the values of a and A. Since 
supj,>Q (5(2;) = Amin{l, a63'^}, condition (A5) is satisfied if 

Q-A 

3-^ <ah< ^A{a) or if {ah < 3"^) A (A < A{a)) . (76) 
A 

This condition reduces to (74) if A = 1; while (76) is never satisfied if A > 1, it holds for a certain 
range of a and 6, if A < 1. For instance, with a = 1/3 we have A{a) = 1/2 and (76) reads b < \/3, 
A < 1/2. 

8 Concluding remarks 

In this paper we have discussed the dynamics of Bianchi type I solutions of the Einstein equations 
with anisotropic matter. The focus of our analysis has been the asymptotic behaviour of solutions 
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in the limit of late times and toward the initial singularity. The matter model was not specified 
explicitly, but only through a set of mild assumptions that are motivated by basic physical consid- 
erations. In this way, our analysis applies to a wide variety of anisotropic matter models including 
matter models as different from each other as collisionless matter, elastic materials and magnetic 
fields. 

A basic assumption on the matter model was to require that the isotropic pressure p and the 
density p obey a linear equation of state p = ivp, where w = const e (—1,1), see Assumption 3. It 
is important to note, however, that this assumption is not necessary. On the contrary, based on 
the results derived in this paper it is relatively simple to treat nonlinear equations of state. Let us 
elaborate. Assume that p is given by a more general function instead of (7), i.e., p = p{n, si, S2, S3). 
In this case we obtain from (5) that w = w{n, si, 82,83) and Wi = Wi{n, si, 82, 83), i = 1,2,3; an 
interesting subcase is p{n, 81, 82, S3) = (p{n)Tp{8i, 82, S3); here, w = w{n) = {d\og(p/dn) — 1 and 
Wi = Wi{n, si,S2, S3) is given by (8), where w is replaced by w(n). The Einstein evolution equations, 
written in terms of the dynamical systems variables of Section 3, decouple into an equation for 
H and a reduced system of equations for the remaining variables, where we can replace x by n, 
since x = n^/^(siS2S3)~^/^. Accordingly, the dynamical system that encodes the dynamics is given 
by (11), which is supplemented by the equation n' = — 3n. When we choose to compactify the 
variable n, i.e., when we replace n by = n/(l + n), the state space of this dynamical system is 
3C X T X (0, 1). If we assume an equation of state such that w{N, si, S2, S3) and Wi{N, si, S2, S3) 
possess well-defined limits as A'' ^ and N ^ 1, we can extend the dynamical system to the 
boundaries Xq = XxT x {0} and A'i=3CxTx{l}of the state space. The dynamical system on 
each of these boundary subsets coincides with the system (11) that we have discussed so extensively 
in this paper. Since the variable A'^ is strictly monotone, the asymptotic dynamics of solutions of 
the dynamical system is associated with the limits A^ ^ and A^ ^ 1. Accordingly, asymptotically, 
the flow of the boundary subsets Aq and Xi (and thus the results of the present paper) constitute 
the key to an understanding of the dynamics of the more general problem with nonlinear equations 
of state. 

It is conceivable that several of our assumptions can be relaxed (however, physics might disapprove), 
which could lead to interesting extensions of the results of this paper. For instance, if Assumption 8 
is removed, then there exist several isotropic states of the matter and thus there might exist Bianchi 
type I solutions that isotropize both toward the past and toward the future. Assumption 1, on the 
other hand, could be replaced by the condition that the energy-momentum tensor depends not only 
on the spatial metric but also on the second fundamental form. In this way it is possible to extend 
the analysis to even larger classes of matter fields. This might lead to valuable generalizations of 
the results presented in this paper. 
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